diff --git a/CMakeLists.txt b/CMakeLists.txt
index b7a9139..2c83353 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -65,7 +65,10 @@ target_compile_definitions(dto PRIVATE DTO_STATS_SUPPORT)
 
 # Link libraries
 set(DTO_STATIC_LIBS "")
-set(DTO_DYNAMIC_LIBS dl)
+# pthread: stats use thread keys and a sampler thread; before glibc 2.34 these
+# are not in libc, so a program that doesn't link libpthread itself would fail
+# to load libdto.
+set(DTO_DYNAMIC_LIBS dl pthread)
 
 # Handle accel-config support
 if(DTO_ACCEL_CONFIG_SUPPORT)
@@ -157,6 +160,11 @@ else()
 endif()
 install(TARGETS dtoctl RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR})
 
+# dto-stats-report summarizes "dtoctl --stats-run"; dtoctl looks for it next to
+# itself first, so copy it beside dtoctl in the build tree too.
+configure_file(dto-stats-report ${CMAKE_CURRENT_BINARY_DIR}/dto-stats-report COPYONLY)
+install(PROGRAMS dto-stats-report DESTINATION ${CMAKE_INSTALL_BINDIR})
+
 # Build and install the optional section 1 manual.
 find_program(PANDOC_EXECUTABLE NAMES pandoc)
 if(PANDOC_EXECUTABLE)
@@ -221,6 +229,33 @@ if(DTO_BUILD_TESTS)
         LABELS "functional"
         TIMEOUT 60
     )
+    # Stats collection (DTO_COLLECT_STATS): a workload of known calls, checked
+    # by a Python script. Runs without DSA (everything is counted as CPU).
+    find_package(Python3 COMPONENTS Interpreter QUIET)
+    add_executable(dto-test-stats-workload tests/test_stats_workload.c)
+    target_compile_options(dto-test-stats-workload PRIVATE -fno-builtin)
+    target_link_libraries(dto-test-stats-workload PRIVATE DTO::dto pthread)
+    if(Python3_Interpreter_FOUND)
+        add_test(NAME stats
+            COMMAND ${Python3_EXECUTABLE} ${CMAKE_CURRENT_SOURCE_DIR}/tests/test_stats.py
+                    $<TARGET_FILE:dto-test-stats-workload> $<TARGET_FILE:dto>)
+        set_tests_properties(stats PROPERTIES
+            LABELS "functional"
+            TIMEOUT 120
+        )
+    else()
+        message(STATUS "stats test: disabled (python3 not found)")
+    endif()
+    # dto-stats-report's arithmetic, on made-up runs
+    if(Python3_Interpreter_FOUND)
+        add_test(NAME stats-report
+            COMMAND ${Python3_EXECUTABLE} ${CMAKE_CURRENT_SOURCE_DIR}/tests/test_stats_report.py
+                    ${CMAKE_CURRENT_BINARY_DIR}/dto-stats-report)
+        set_tests_properties(stats-report PROPERTIES
+            LABELS "functional"
+            TIMEOUT 60
+        )
+    endif()
 
     # ---- A/B performance regression test ----
     # The baseline object is compiled from unguarded upstream sources that
diff --git a/Makefile b/Makefile
index 25b649a..15467f5 100644
--- a/Makefile
+++ b/Makefile
@@ -16,7 +16,7 @@ MANDIR ?= /usr/share/man
 DTO_VERSION := $(shell grep -o 'DTO VERSION [0-9.]*' CMakeLists.txt | cut -d' ' -f3)
 
 libdto: dto.c
-	gcc -shared -fPIC -Wl,-soname,libdto.so dto.c $(DML_LIB_CXX) -DDTO_STATS_SUPPORT -DDTO_ACCEL_CONFIG_SUPPORT -DDTO_NUMA_SUPPORT -o libdto.so.1.0 -laccel-config -ldl -lnuma -mwaitpkg
+	gcc -shared -fPIC -Wl,-soname,libdto.so dto.c $(DML_LIB_CXX) -DDTO_STATS_SUPPORT -DDTO_ACCEL_CONFIG_SUPPORT -DDTO_NUMA_SUPPORT -o libdto.so.1.0 -laccel-config -ldl -lnuma -lpthread -mwaitpkg
 
 libdto_nostats: dto.c
 	gcc -shared -fPIC -Wl,-soname,libdto.so dto.c $(DML_LIB_CXX) -DDTO_ACCEL_CONFIG_SUPPORT -DDTO_NUMA_SUPPORT -o libdto.so.1.0 -laccel-config -ldl -lnuma -mwaitpkg
@@ -43,6 +43,7 @@ install: dtoctl install-man
 	ln -sf /usr/lib64/libdto.so.1.0 /usr/lib64/libdto.so.1
 	ln -sf /usr/lib64/libdto.so.1.0 /usr/lib64/libdto.so
 	cp dtoctl /usr/bin/
+	cp dto-stats-report /usr/bin/
 
 install-man: man
 	@if [ -f doc/dtoctl.1 ]; then \
diff --git a/README.md b/README.md
index 71e1424..bed25a0 100644
--- a/README.md
+++ b/README.md
@@ -50,9 +50,12 @@ dto-4-dsa.conf:  An example json config file for configuring DSAs
 Following environment variables control the behavior of DTO library:
 	DTO_USESTDC_CALLS=0/1, 1 (uses std c memory functions only), 0 (uses DSA along with std c lib call; in case of DSA page fault - reverts to std c lib call). Default is 0.
         DTO_COLLECT_STATS=0/1, 1 (enables stats collection - #of operations, avg latency for each API, etc.>, 0 (disables stats collection).
-                                Should be enabled for debugging/profiling only, not for perf evaluation (enabling it slows down the workload). Default is 0.
-        DTO_STATS_FILE=<path> when DTO_COLLECT_STATS=1, redirect histogram output to the provided regular file instead of standard output.
-                                The file is truncated when DTO loads; failures fall back to standard output and emit an error.
+                                Adds about 40 ns per call (two clock reads), whatever the number of threads. Default is 0. See Statistics below.
+        DTO_STATS_FILE=<path> when DTO_COLLECT_STATS=1, write stats to the provided regular file instead of standard output. Each %p in
+                                the path is replaced by the process ID. Records are appended, not truncated: processes started from this one
+                                append to the same file unless the path has a %p. Failures fall back to standard output and emit an error.
+        DTO_STATS_FORMAT=text/json  text (default): the tables shown below; json: JSON Lines, one record per line (see Statistics).
+        DTO_STATS_INTERVAL_MS=N  when N > 0, also write what happened in every N ms, from a background thread. Default is 0.
 	DTO_WAIT_METHOD=<yield,busypoll,umwait,tpause> (specifies the method to use while waiting for DSA to complete operation, default is busypoll)
 	DTO_MIN_BYTES=xxxx (specifies minimum size of API call needed for DSA operation execution, default is 65536 bytes)
 	DTO_CPU_SIZE_FRACTION=0.xx (specifies fraction of job performed by CPU, in parallel to DSA). Default is 0.00
@@ -192,6 +195,89 @@ rather than `dtoctl sudo ...`, which loses `LD_PRELOAD` and the `DTO_*`
 variables, and pass settings as options: `sudo` does not pass on variables
 exported in your shell.
 
+### Measuring the speedup from DSA
+
+`dtoctl --stats-run` runs a program CPU-only (`DTO_USESTDC_CALLS=1`) and then
+with DSA, collecting DTO's JSON stats (see [Statistics](#statistics)) from each
+run into a directory, and `dto-stats-report` (installed with `dtoctl`) compares
+them. With `--repeat=N` each way runs N times, alternating, and the report uses
+the median wall time and, for the CPU and mem* times, the run with the median
+CPU time; `--perf` adds one more run each way under `perf record`.
+
+```bash
+sudo taskset -c 2 dtoctl -w busypoll --no-auto-adjust -c 0 --stats-run=app-stats --repeat=3 --perf ./app
+```
+
+```
+DTO stats run: app-stats
+  command: ./app
+  settings: DTO_AUTO_ADJUST_KNOBS=0 DTO_CPU_SIZE_FRACTION=0 DTO_WAIT_METHOD=busypoll
+  timed runs: 3 CPU-only, 3 with DSA (medians)
+  DTO in the DSA run: 1 WQs, wait busypoll, min_bytes 65536, cpu_size_fraction 0.0
+
+                                     CPU only       with DSA   speedup
+  wall time                            1.66 s         1.38 s     1.20x
+  CPU time (user+sys)                  1.65 s         1.38 s     1.20x
+    in mem* calls                   782.69 ms      493.11 ms     1.59x
+    outside mem* calls              868.09 ms      885.71 ms     0.98x
+  wall time, min-max              1.66-1.67 s    1.38-1.40 s
+
+  Amdahl's law
+    f, share of CPU time in mem* calls (CPU-only)     47.4%
+    s, mem* speedup with DSA (time per byte)          1.59x
+    predicted speedup, 1 / ((1 - f) + f / s)          1.21x
+    measured speedup, wall time                       1.20x
+    limit if mem* took no time, 1 / (1 - f)           1.90x
+
+  per operation       calls       bytes   avg size   CPU only   with DSA  speedup  on DSA       f
+  memcpy             20.00k     4.9 GiB  256.0 KiB  596.44 ms  180.24 ms    3.31x   99.4%   36.1%
+  memset             20.00k     4.9 GiB  256.0 KiB   92.17 ms  100.59 ms    0.92x  100.0%    5.6%
+  memcmp             20.00k     4.9 GiB  256.0 KiB   94.09 ms  212.28 ms    0.44x  100.0%    5.7%
+    (on DSA: share of the bytes in calls sent to DSA, of which the CPU does cpu_size_fraction;
+     f: share of the CPU-only run's CPU time)
+    calls DSA left to the CPU, per run: page fault 128, retry 0, other 0
+
+  by size (most time)       calls   CPU only   with DSA  speedup  on DSA
+  256K-320K                60.00k   13.04 us    8.22 us    1.59x   99.8%
+    (average time per call; on DSA: share of calls sent to DSA)
+
+  profile (perf, % of CPU time)        CPU only   with DSA
+    memcpy/memmove                        35.0%       0.3%
+    memset                                 6.0%       0.0%
+    memcmp                                 5.3%       0.0%
+    CRC                                    9.8%      12.5%
+    DTO library                            0.3%      28.7%
+    other user code                       31.2%      43.8%
+    other kernel code                     12.4%      14.8%
+
+  note: 10% of the CPU-only run is CRC code: DSA can compute CRC32C, but DTO doesn't intercept CRC functions, so offloading it needs an explicit call
+```
+
+Here a single thread, on a core of the DSA's NUMA node, copies 256 KiB at a
+time out of a 1 GiB buffer, so the data is not in the CPU caches, compares the
+copy with the source, clears it, computes a CRC32C in software, and does other
+work. The CPU-only run spent f = 47.4% of its CPU time in mem\* calls, and
+those calls took s = 1.59x less time per byte with DSA, so Amdahl's law
+predicts 1 / ((1 - f) + f / s) = 1.21x for the whole program; it measured
+1.20x. At best, mem\* calls taking no time at all, it could run 1.90x faster.
+The copies were 3.31x faster with DSA, but the compares, which read the copy
+right after it is made, were slower: a copy on the CPU leaves the data in its
+caches and a copy by DSA does not. The auto-tuner is off in this run
+(`--no-auto-adjust -c 0`); with it on, it gave part of each copy to the CPU and
+the program ran only 1.02x faster.
+
+The report also shows each operation and the sizes with the most time, so you
+can see where DSA wins and set `DTO_MIN_BYTES` accordingly, and notes anything
+that makes the comparison doubtful. Among them: no call reaching DSA (no usable
+work queue: run with `sudo`), and runs that did different amounts of work, as
+a benchmark that runs for a fixed time does (speedups are then per byte). The
+perf profile splits each run's CPU time into copies, memset, memcmp, CRC,
+kernel copies, DTO itself (here, busy-polling for completions) and the rest. It
+also counts mem\* calls DTO cannot intercept, such as those made inside libc.
+Run `dto-stats-report app-stats` again at any time, or with `--json` for the
+numbers. Amdahl's estimate assumes the program keeps its CPUs busy; one that
+waits for I/O gains less in wall time.
+
 Run `dtoctl --help` for the option list, or `man 1 dtoctl` for the full manual,
 which documents each option and gives a ready-to-run command for each of the
 usage models above.
@@ -310,63 +396,76 @@ to the DSA device.
 
 	<CBENCH_DIR>/cachebench -json_test_config <json file> --progress_stats_file=dto.log --report_api_latency
 
-4. Sample histogram given below (generated using DTO_COLLECT_STATS=1)
+4. Sample stats output (DTO_COLLECT_STATS=1; see [Statistics](#statistics))
 	i. Numbers under columns set, cpy, mov, and cmp show number of API calls or per-API completion latency for memset, memcpy, memmove, and memcmp respectively.
 	ii. Numbers in bytes column show total bytes processed across all 4 API calls
 
 ------------------------------------------------------------------------------------------------------------------------------------------
 
-DTO Run Time: 22826 ms
+```
+DTO Run Time: 8 ms
 
 ******** Number of Memory Operations ********
-                     <*************** stdc calls    ***************> <*************** dsa (success) ***************> <*************** dsa (failed)  ***************> <***** failure reason *****>
-Byte Range        -- set      cpy      mov      cmp      bytes        set      cpy      mov      cmp      bytes        set      cpy      mov      cmp      bytes        Retries PFs    Others
-       0-4095     -- 24616626 12319310 7676731  244777209 4531255821  0        0        0        0        0            0        0        0        0        0            0      0      0
-    4096-8191     -- 3        290      14       0        1842507      0        0        0        0        0            0        0        0        0        0            0      0      0
-    8192-12287    -- 0        2        6        0        45579        64       176      2        0        2545968      0        2        6        0        25276        0      8      0
-   12288-16383    -- 1        1        0        0        20756        1        137      0        0        1965343      1        1        0        0        11172        0      2      0
-   16384-20479    -- 0        1        3        0        41974        0        126      2        0        2328434      0        1        3        0        24337        0      4      0
-   32768-36863    -- 0        0        0        0        0            0        67       2        0        2378165      0        0        0        0        0            0      0      0
-   53248-57343    -- 0        0        0        0        0            0        46       0        0        2545890      0        0        0        0        0            0      0      0
-   98304-102399   -- 0        0        0        0        0            0        13       0        0        1309627      0        0        0        0        0            0      0      0
-  106496-110591   -- 0        1        0        0        31637        0        8        0        0        869880       0        1        0        0        77586        0      1      0
-  118784-122879   -- 0        0        0        0        0            0        9        0        0        1083519      0        0        0        0        0            0      0      0
-  131072-135167   -- 0        0        2        0        175638       0        5        0        0        665475       0        0        2        0        86506        0      2      0
-  143360-147455   -- 0        0        0        0        0            0        1        0        0        143399       0        0        0        0        0            0      0      0
-  147456-151551   -- 0        0        0        0        0            0        1        0        0        148647       0        0        0        0        0            0      0      0
-  159744-163839   -- 0        0        0        0        0            0        6        0        0        970186       0        0        0        0        0            0      0      0
-  212992-217087   -- 0        0        11       0        1580248      0        4        2        0        1292252      0        0        11       0        779384       0      11     0
-  278528-282623   -- 0        0        0        0        0            0        1        0        0        279719       0        0        0        0        0            0      0      0
-  466944-471039   -- 0        0        0        0        0            0        1        0        0        467783       0        0        0        0        0            0      0      0
-  528384-532479   -- 0        0        0        0        0            0        1        0        0        529927       0        0        0        0        0            0      0      0
-  720896-724991   -- 0        0        0        0        0            1        0        0        0        720896       0        0        0        0        0            0      0      0
-  999424-1003519  -- 0        0        0        0        0            0        1        0        0        1002439      0        0        0        0        0            0      0      0
-   >=2093056      -- 0        1        0        0        1975911      0        0        0        0        0            0        1        0        0        973209       0      1      0
+                           <*************** stdc calls    ***************> <*************** dsa (success) ***************> <*************** dsa (failed)  ***************> <***** failure reason *****>
+Byte Range              -- set      cpy      mov      cmp      bytes        set      cpy      mov      cmp      bytes        set      cpy      mov      cmp      bytes        Retry  PFs    Others
+      65536-81919       -- 0        0        0        0        0            200      0        0        0        13926400     0        0        0        0        0            0      0      0
+      98304-114687      -- 0        0        200      0        19662200     0        0        0        0        0            0        0        0        0        0            0      0      0
+     131072-163839      -- 0        0        0        0        0            0        200      0        0        26239000     0        0        0        0        0            0      0      0
+     262144-327679      -- 0        0        0        0        0            0        0        0        200      52441600     0        0        0        0        0            0      0      0
 
 ******** Average Memory Operation Latency (us)  ********
-                     <******** stdc calls    ********> <******** dsa (success) ********> <******** dsa (failed)  ********> 
-Byte Range        -- set      cpy      mov      cmp      set      cpy      mov      cmp      set      cpy      mov      cmp
-       0-4095     -- 0.01     0.02     0.01     0.04     0        0        0        0        0        0        0        0
-    4096-8191     -- 0.07     0.42     0.47     0        0        0        0        0        0        0        0        0
-    8192-12287    -- 0        1.29     0.47     0        0.84     1.19     1.82     0        0        56.18    2.54     0
-   12288-16383    -- 2.90     3.34     0        0        7.49     1.26     0        0        2.07     1.93     0        0
-   16384-20479    -- 0        0.27     3.24     0        0        1.37     1.86     0        0        84.92    2.47     0
-   32768-36863    -- 0        0        0        0        0        1.87     2.91     0        0        0        0        0
-   53248-57343    -- 0        0        0        0        0        2.45     0        0        0        0        0        0
-   98304-102399   -- 0        0        0        0        0        3.76     0        0        0        0        0        0
-  106496-110591   -- 0        229.18   0        0        0        4.19     0        0        0        3.65     0        0
-  118784-122879   -- 0        0        0        0        0        5.18     0        0        0        0        0        0
-  131072-135167   -- 0        0        25.00    0        0        5.92     0        0        0        0        14.39    0
-  143360-147455   -- 0        0        0        0        0        5.11     0        0        0        0        0        0
-  147456-151551   -- 0        0        0        0        0        6.16     0        0        0        0        0        0
-  159744-163839   -- 0        0        0        0        0        7.65     0        0        0        0        0        0
-  212992-217087   -- 0        0        36.10    0        0        8.15     6.59     0        0        0        18.07    0
-  278528-282623   -- 0        0        0        0        0        15.47    0        0        0        0        0        0
-  466944-471039   -- 0        0        0        0        0        15.91    0        0        0        0        0        0
-  528384-532479   -- 0        0        0        0        0        18.65    0        0        0        0        0        0
-  720896-724991   -- 0        0        0        0        17.71    0        0        0        0        0        0        0
-  999424-1003519  -- 0        0        0        0        0        28.42    0        0        0        0        0        0
-   >=2093056      -- 0        422.15   0        0        0        0        0        0        0        272.87   0        0
+                           <******** stdc calls    ********> <******** dsa (success) ********> <******** dsa (failed)  ********>
+Byte Range              -- set      cpy      mov      cmp      set      cpy      mov      cmp      set      cpy      mov      cmp
+      65536-81919       -- 0        0        0        0        8.06     0        0        0        0        0        0        0
+      98304-114687      -- 0        0        4.18     0        0        0        0        0        0        0        0        0
+     131072-163839      -- 0        0        0        0        0        17.02    0        0        0        0        0        0
+     262144-327679      -- 0        0        0        0        0        0        0        34.50    0        0        0        0
+```
+
+## Statistics
+
+With `DTO_COLLECT_STATS=1`, DTO times and counts every memset, memcpy,
+memmove and memcmp it intercepts, by operation, by path (CPU, DSA, or DSA
+failed and finished on the CPU) and by size, and writes a summary when the
+process exits (or on Ctrl-C). Each thread counts on its own, so collection
+costs about 40 ns per call however many threads there are, and the counts
+are summed for output.
+
+The size histogram has a bucket for calls under 64 bytes, then four buckets
+per power of two up to 4 GiB (65536-81919, 81920-98303, ...), and one for
+anything larger. A call that DSA completed only partly (for example on a
+page fault) is counted once as `dsa (failed)` and once as a CPU call for the
+remainder.
+
+With `DTO_STATS_INTERVAL_MS=N`, a background thread also writes, every N ms,
+what happened since the previous record, so activity can be followed over
+time. The interval records add up to the summary's totals.
+
+With `DTO_STATS_FORMAT=json` the output is JSON Lines, one record per line:
+
+```
+{"type":"interval","pid":183654,"t_ms":14.333,"dt_ms":14.333,"ops":{"memcpy":{"dsa":{"calls":8,"bytes":1049560,"ns":96273}},"memmove":{"cpu":{"calls":8,"bytes":786488,"ns":31601}}},"failures":{"retry":0,"page_fault":0,"other":0}}
+{"type":"summary","pid":183654,"program":"prog","run_ms":23.395,
+ "config":{"stdc_only":0,"wait_method":"busypoll","min_bytes":65536,"cpu_size_fraction":0.00,"auto_adjust":1,"numa_aware":0,"wqs":1},
+ "totals":{"memcpy":{"dsa":{"calls":12,"bytes":1574340,"ns":138591}},...},
+ "failures":{"retry":0,"page_fault":0,"other":0},
+ "histogram":[{"min":131072,"max":163839,"ops":{"memcpy":{"dsa":{"calls":12,"bytes":1574340,"ns":138591}}},"failures":{...}},...]}
+```
+
+(The summary is a single line; it is wrapped here.) `ops` maps `memset`,
+`memcpy`, `memmove` and `memcmp` to the paths they took, `cpu`, `dsa` and
+`dsa_failed`, each with the number of calls, bytes and total nanoseconds;
+operations and paths with no calls are left out. `t_ms` and `run_ms` count
+from DTO's start, and `config` holds the settings at exit, including the
+auto-tuned `min_bytes` and `cpu_size_fraction`.
+
+To compare a run with DSA against a CPU-only run, and see how close the speedup
+comes to what Amdahl's law predicts, use `dtoctl --stats-run` (see
+[Measuring the speedup from DSA](#measuring-the-speedup-from-dsa)).
+
+Every process that loads DTO with stats enabled writes its own records,
+tagged with its `pid`. They share `DTO_STATS_FILE` (records are appended
+whole) unless the path contains `%p`, which gives each process its own file.
 
 ## Usage Notes
 
diff --git a/TESTING.md b/TESTING.md
index db6b543..bf2ef7b 100644
--- a/TESTING.md
+++ b/TESTING.md
@@ -86,8 +86,10 @@ queues DTO uses share the same setting.
 Checks for the `dtoctl` wrapper: every option sets its `DTO_*` variable,
 invalid values are rejected, the library path is resolved and validated,
 `LD_PRELOAD` and inherited variables are handled, and exit statuses follow
-`env(1)`. The last check runs a program under `dtoctl` and reads libdto's own
-configuration trace. No DSA hardware is needed.
+`env(1)`. `--stats-run` is checked from its `--dry-run` plan and by running a
+small shell command: each run's environment and stats file, `runs.jsonl`, the
+report, and stopping at a failed run. The last check runs a program under
+`dtoctl` and reads libdto's own configuration trace. No DSA hardware is needed.
 
 ```bash
 ctest -R dtoctl
@@ -95,6 +97,25 @@ ctest -R dtoctl
 tests/test_dtoctl.sh build/dtoctl build/libdto.so
 ```
 
+## Stats Tests
+
+`stats` runs `dto-test-stats-workload`, which makes a known number of calls
+of known sizes, with `DTO_COLLECT_STATS` set, and checks the output: totals
+and histogram buckets, interval records adding up to the totals, calls from
+200 short-lived threads all counted without their stats memory piling up,
+separate records for a forked child, `%p` file names, and the text format.
+It needs Python 3 and no DSA hardware.
+
+```bash
+ctest -R stats
+```
+
+`stats-report` checks `dto-stats-report` on made-up runs with known answers
+(Amdahl's law, per-operation and per-size speedups, medians, the notes it
+makes), sorts made-up perf samples into categories, and checks that it finds
+memcpy, memset and memcmp by address in this machine's libc. It needs Python 3,
+and neither DSA nor perf.
+
 ## Performance Tests
 
 The performance test (`perf`) detects regressions in libdto by
diff --git a/doc/dtoctl.md b/doc/dtoctl.md
index afca4ed..be73825 100644
--- a/doc/dtoctl.md
+++ b/doc/dtoctl.md
@@ -13,6 +13,8 @@ dtoctl [_options_] COMMAND [_command-options_]
 
 dtoctl \--dry-run [_options_] [COMMAND [_command-options_]]
 
+dtoctl \--stats-run[=_DIR_] [\--repeat=_N_] [\--perf[=_PERF_]] [_options_] COMMAND [_command-options_]
+
 
 DESCRIPTION
 ===========
@@ -160,13 +162,24 @@ Other toggles
 Statistics and logging
 ----------------------
 \--stats
-:   Enable stats collection.  For debugging/profiling only, as it slows down the
-    workload.  Sets **DTO_COLLECT_STATS=1**.  Default: off.
+:   Enable stats collection: calls, bytes and time for each operation, by
+    path (CPU, DSA, DSA failed) and by size.  Costs about 40 ns per call.
+    Sets **DTO_COLLECT_STATS=1**.  Default: off.
 
 \--stats-file=_PATH_
-:   Write the stats histogram to _PATH_ instead of standard output.  Implies
-    `--stats`, since DTO ignores the file otherwise.  Sets **DTO_STATS_FILE**
-    and **DTO_COLLECT_STATS=1**.  Default: standard output.
+:   Write the stats to _PATH_ instead of standard output; each `%p` in _PATH_
+    becomes the process ID.  The file is appended to.  Implies `--stats`,
+    since DTO ignores the file otherwise.  Sets **DTO_STATS_FILE** and
+    **DTO_COLLECT_STATS=1**.  Default: standard output.
+
+\--stats-format=_FORMAT_
+:   _text_ for tables, or _json_ for JSON Lines (one record per line).
+    Implies `--stats`.  Sets **DTO_STATS_FORMAT**.  Default: _text_.
+
+\--stats-interval=_MS_
+:   Also write the stats for every _MS_ milliseconds, from a background
+    thread, as well as the summary at exit.  Implies `--stats`.  Sets
+    **DTO_STATS_INTERVAL_MS**.  Default: summary only.
 
 \--log-file=_PATH_
 :   Redirect DTO output to _PATH_._program_._pid_, one file per process.
@@ -176,6 +189,28 @@ Statistics and logging
 :   Verbosity of DTO logging: _0_, _1_, or _2_.  Sets **DTO_LOG_LEVEL**.
     Default: _0_.
 
+Stats run
+---------
+\--stats-run[=_DIR_]
+:   Measure how much DSA speeds `COMMAND` up.  Run it CPU-only
+    (**DTO_USESTDC_CALLS=1**), then with DSA (**DTO_USESTDC_CALLS=0**),
+    collecting JSON stats from each run into _DIR_, then summarize them with
+    **dto-stats-report** (see STATS RUN).  _DIR_ must be new or empty; the
+    default is a new `dto-stats-<date>-<time>` in the current directory.
+    Give _DIR_ with `=`: `--stats-run=out`.  The other options apply to every
+    run; `--stdc-only`, `--stats-file` and `--stats-format` cannot be used,
+    since each run sets them.
+
+\--repeat=_N_
+:   With `--stats-run`, run each way _N_ times, alternating; the report takes
+    the median wall time, and the CPU and mem* times of the run with the
+    median CPU time.  Default: 1.
+
+\--perf[=_PERF_]
+:   With `--stats-run`, also run each way once more under
+    `PERF record -e cpu-clock`, to show where the CPU time went.  _PERF_ is
+    the **perf** to run.  Default: **perf** from **PATH**.
+
 -?, \--help
 :   Print help message and list of options with description.
 
@@ -205,6 +240,51 @@ status is that of `COMMAND`.  Otherwise:
 
 `--dry-run`, `--help`, and `--version` exit with 0.
 
+With `--stats-run`, a run that fails stops the rest, and **dtoctl** exits with
+its status (128 + the signal number if it was killed); otherwise the exit
+status is that of **dto-stats-report**.
+
+
+STATS RUN
+=========
+`--stats-run` answers two questions: how much faster is `COMMAND` with DSA,
+and how much faster could it be?
+
+Each run's DTO stats give the time spent in the intercepted calls, and
+**dtoctl** records each run's wall time and CPU time (user + system, for
+`COMMAND` and all its children).  **dto-stats-report** then works out:
+
+f
+:   the share of the CPU-only run's CPU time spent in mem\* calls;
+
+s
+:   how much faster the mem\* calls were with DSA, per byte;
+
+predicted speedup
+:   by Amdahl's law, 1 / ((1 - f) + f / s), set against the measured
+    speedup in wall time;
+
+limit
+:   1 / (1 - f), the speedup if the mem\* calls took no time at all.
+
+It also lists each operation and the sizes with the most time, CPU-only and
+with DSA, and notes where DSA lost, where the runs differed, and why no call
+may have reached DSA.  With `--perf` it divides each run's CPU time into
+copies, memset, memcmp, CRC, kernel copies, DTO itself and other code, from
+the profile, which also catches mem\* calls DTO cannot intercept.
+
+Amdahl's estimate assumes `COMMAND` is busy on the CPU; one that waits for
+I/O or for other threads gains less in wall time.  A benchmark that runs for
+a fixed time does more work when faster: the report says so when the runs
+copied different amounts, and compares speed per byte.
+
+The directory holds `runs.jsonl` (the command, settings, and each run's times
+and exit status), `<run>.jsonl` (DTO's stats for each run, from every process
+it started), and with `--perf`, `<mode>-perf.data`.  Run
+**dto-stats-report** _DIR_ again at any time, with `--json` for the numbers
+as JSON.  DTO writes the stats files readable by their owner only, so a run
+made with **sudo** needs **sudo** to report on.
+
 
 NOTES
 =====
@@ -255,6 +335,10 @@ stats into a file:
     $ dtoctl -l ./libdto.so.1.0 -q "wq0.0;wq2.0;wq4.0;wq6.0" \
           --stats-file=/tmp/dto-stats.txt ./prog
 
+Measure the speedup from DSA, three runs each way, then profile:
+
+    $ sudo dtoctl --stats-run=prog-stats --repeat=3 --perf ./prog
+
 Check what would run, ignoring anything exported in the shell:
 
     $ dtoctl --dry-run -i -w umwait -b 16K ./prog
@@ -267,4 +351,4 @@ Check what would run, ignoring anything exported in the shell:
 
 SEE ALSO
 ========
-**ld.so**(8), **accel-config**(1)
+**ld.so**(8), **accel-config**(1), **perf-record**(1)
diff --git a/dto-stats-report b/dto-stats-report
new file mode 100755
index 0000000..75b51d4
--- /dev/null
+++ b/dto-stats-report
@@ -0,0 +1,685 @@
+#!/usr/bin/env python3
+# SPDX-License-Identifier: MIT
+"""Summarize a "dtoctl --stats-run": how much faster the program ran with DSA
+than CPU-only, and how that compares with Amdahl's law, using the time the
+CPU-only run spent in mem* calls and how much faster those calls were with
+DSA. With a perf profile of each run, also shows where the CPU time went:
+copies, memset, memcmp, CRC, kernel copies, DTO itself.
+
+usage: dto-stats-report [--json] [--perf PERF] [--top N] DIR
+"""
+
+import argparse
+import bisect
+import ctypes
+import json
+import os
+import re
+import statistics
+import struct
+import subprocess
+import sys
+
+OPS = ("memcpy", "memmove", "memset", "memcmp")
+PATHS = ("cpu", "dsa", "dsa_failed")
+FAILS = ("retry", "page_fault", "other")
+FIELDS = ("calls", "bytes", "ns")
+MODES = ("cpu", "dsa")
+
+notes = []
+
+
+def note(msg):
+    notes.append(msg)
+
+
+# ---- loading ---------------------------------------------------------------
+
+def load_jsonl(path):
+    recs = []
+    try:
+        with open(path, encoding="utf-8", errors="replace") as f:
+            for n, line in enumerate(f, 1):
+                if not line.strip():
+                    continue
+                try:
+                    recs.append(json.loads(line))
+                except json.JSONDecodeError:
+                    note(f"{os.path.basename(path)} line {n} is not valid JSON "
+                         "(cut short?); skipped")
+    except FileNotFoundError:
+        pass
+    except PermissionError:
+        # DTO writes stats files readable by their owner only
+        sys.exit(f"dto-stats-report: cannot read {path}: permission denied; the "
+                 "files belong to the user that ran dtoctl (run this with sudo?)")
+    return recs
+
+
+def zero_ops():
+    return {op: {p: dict.fromkeys(FIELDS, 0) for p in PATHS} for op in OPS}
+
+
+def add_ops(acc, ops, scale=1.0):
+    for op, paths in ops.items():
+        if op not in acc:
+            continue
+        for p, v in paths.items():
+            if p in acc[op]:
+                for k in FIELDS:
+                    acc[op][p][k] += v.get(k, 0) * scale
+
+
+class RunStats:
+    """DTO's stats for one or more runs of one mode, summed over every
+    process each run started, averaged over the runs."""
+
+    def __init__(self):
+        self.ops = zero_ops()
+        self.fails = dict.fromkeys(FAILS, 0)
+        self.hist = {}          # (min, max) -> ops
+        self.configs = []
+        self.processes = 0
+
+    def add_run(self, summaries, scale):
+        for s in summaries:
+            add_ops(self.ops, s.get("totals", {}), scale)
+            for k in FAILS:
+                self.fails[k] += s.get("failures", {}).get(k, 0) * scale
+            for b in s.get("histogram", []):
+                key = (b["min"], b.get("max"))
+                add_ops(self.hist.setdefault(key, zero_ops()), b.get("ops", {}), scale)
+            self.configs.append(s.get("config", {}))
+            self.processes += scale
+
+    def op_total(self, op, field, paths=PATHS):
+        return sum(self.ops[op][p][field] for p in paths)
+
+    def total(self, field, paths=PATHS):
+        return sum(self.op_total(op, field, paths) for op in OPS)
+
+
+def load_dir(d):
+    recs = load_jsonl(os.path.join(d, "runs.jsonl"))
+    if not recs or recs[0].get("type") != "stats_run":
+        sys.exit(f"dto-stats-report: {d} has no runs.jsonl from dtoctl --stats-run")
+    header, runs = recs[0], [r for r in recs[1:] if r.get("type") == "run"]
+    no_file, no_summary = [], []
+    for r in runs:
+        if r.get("exit_status", 0) or "signal" in r:
+            how = (f"was killed by signal {r['signal']}" if "signal" in r
+                   else f"exited with status {r['exit_status']}")
+            if r.get("profiled"):
+                note(f"profiled run {r['name']} {how}; the timed runs are "
+                     "reported without a profile")
+            else:
+                note(f"run {r['name']} {how}; later runs did not happen")
+        path = os.path.join(d, r["stats"])
+        r["summaries"] = [s for s in load_jsonl(path) if s.get("type") == "summary"]
+        if not r["summaries"] and r.get("exit_status") == 0:
+            (no_summary if os.path.exists(path) else no_file).append(r["name"])
+    if no_file:
+        note(f"no DTO stats from {runs_list(no_file)}: is {header.get('library')} "
+             "new enough to write JSON stats, and is the program dynamically linked "
+             "and not setuid?")
+    if no_summary:
+        note(f"DTO stats from {runs_list(no_summary)} have no summary: did the "
+             "program end without running exit handlers (_exit, or a signal)?")
+    return header, runs
+
+
+def runs_list(names):
+    return ("run " if len(names) == 1 else "runs ") + ", ".join(names)
+
+
+# ---- Amdahl ----------------------------------------------------------------
+
+def mode_summary(runs):
+    timed = [r for r in runs if r.get("exit_status") == 0 and r["summaries"]]
+    if not timed:
+        return None
+    # The CPU time, the mem* time and the histogram all come from one run,
+    # the one with the median CPU time, so that f and "outside mem* calls"
+    # compare like with like: a mean mem* time against a median CPU time
+    # lets one preempted run push f past 100%. The wall time is the median.
+    by_cpu = sorted(timed, key=lambda r: r["user_s"] + r["sys_s"])
+    rep = by_cpu[(len(by_cpu) - 1) // 2]
+    st = RunStats()
+    st.add_run(rep["summaries"], 1.0)
+    walls = [r["wall_s"] for r in timed]
+    return {
+        "runs": len(timed),
+        "wall_s": statistics.median(walls),
+        "wall_range": (min(walls), max(walls)),
+        "cpu_s": rep["user_s"] + rep["sys_s"],
+        "mem_s": st.total("ns") / 1e9,
+        "bytes": st.total("bytes"),
+        "stats": st,
+    }
+
+
+def ratio(a, b):
+    return a / b if b else None
+
+
+def per_byte_speedup(ns_cpu, bytes_cpu, ns_dsa, bytes_dsa):
+    """Speedup of time per byte, so that runs that did different amounts
+    of work (a benchmark that runs for a fixed time) still compare."""
+    if not (ns_cpu and bytes_cpu and ns_dsa and bytes_dsa):
+        return None
+    return (ns_cpu / bytes_cpu) / (ns_dsa / bytes_dsa)
+
+
+def analyze(header, runs, top):
+    res = {"command": header.get("command", []), "library": header.get("library"),
+           "host": header.get("host"), "env": header.get("env", {})}
+    m = {mode: mode_summary([r for r in runs if r["mode"] == mode
+                             and not r.get("profiled")]) for mode in MODES}
+    res["modes"] = {k: (None if v is None else {x: v[x] for x in v if x != "stats"})
+                    for k, v in m.items()}
+    c, d = m["cpu"], m["dsa"]
+    if c is None or d is None:
+        note("need a successful CPU-only run and a successful DSA run, with stats, "
+             "to compare")
+        return res, m
+    cs, ds = c["stats"], d["stats"]
+
+    # how the DSA run used DSA
+    cfg = ds.configs[-1] if ds.configs else {}
+    res["dsa_config"] = cfg
+    if ds.total("calls", ("dsa", "dsa_failed")) == 0:
+        if any(x.get("stdc_only") for x in ds.configs):
+            note("no call went to DSA in the DSA run: DTO found no usable work queue "
+                 "and used the CPU (are the WQs enabled, and is the program allowed to "
+                 "open them? try sudo)")
+        else:
+            note(f"no call went to DSA in the DSA run: every call was smaller than "
+                 f"min_bytes ({cfg.get('min_bytes', '?')}) or its operation was not "
+                 "offloaded")
+    if c["bytes"] and d["bytes"] and abs(d["bytes"] / c["bytes"] - 1) > 0.02:
+        note(f"the runs did different amounts of mem* work ({fmt_bytes(c['bytes'])} "
+             f"CPU-only, {fmt_bytes(d['bytes'])} with DSA), so the wall-time speedup "
+             "compares unlike runs; if the program runs for a fixed time, compare "
+             f"mem* throughput: {fmt_bytes(c['bytes'] / c['wall_s'])}/s CPU-only, "
+             f"{fmt_bytes(d['bytes'] / d['wall_s'])}/s with DSA")
+
+    f = ratio(c["mem_s"], c["cpu_s"])
+    if f is not None and f > 1:
+        note(f"the CPU-only run's mem* calls took {f:.0%} of its CPU time: threads "
+             "were preempted or descheduled during calls, so f is taken as 100%")
+        f = 1.0
+    s = per_byte_speedup(c["mem_s"], c["bytes"], d["mem_s"], d["bytes"])
+    amdahl = {"f": f, "s": s, "predicted": None,
+              "limit": None if f is None else (1 / (1 - f) if f < 1 else None)}
+    if f is not None and s:
+        amdahl["predicted"] = 1 / ((1 - f) + f / s)
+    res["amdahl"] = amdahl
+    res["measured"] = {"wall": ratio(c["wall_s"], d["wall_s"]),
+                       "cpu": ratio(c["cpu_s"], d["cpu_s"])}
+    res["outside_mem_s"] = {"cpu": c["cpu_s"] - c["mem_s"],
+                            "dsa": d["cpu_s"] - d["mem_s"]}
+
+    ops = {}
+    for op in OPS:
+        calls, nbytes = cs.op_total(op, "calls"), cs.op_total(op, "bytes")
+        failed = cs.op_total(op, "calls", ("dsa_failed",))
+        if not (calls or ds.op_total(op, "calls")):
+            continue
+        dbytes = ds.op_total(op, "bytes")
+        ops[op] = {
+            # a call DSA left partly done is counted twice: dsa_failed, then cpu
+            "calls": calls - failed,
+            "bytes": nbytes,
+            "cpu_only_s": cs.op_total(op, "ns") / 1e9,
+            "with_dsa_s": ds.op_total(op, "ns") / 1e9,
+            "speedup": per_byte_speedup(cs.op_total(op, "ns"), nbytes,
+                                        ds.op_total(op, "ns"), dbytes),
+            "on_dsa": ratio(ds.op_total(op, "bytes", ("dsa",)), dbytes),
+            "share": ratio(cs.op_total(op, "ns") / 1e9, c["cpu_s"]),
+            "dsa_paths": {p: ds.ops[op][p] for p in PATHS},
+        }
+    res["ops"] = ops
+    res["dsa_failures"] = ds.fails
+
+    buckets = []
+    for key in sorted(set(cs.hist) | set(ds.hist)):
+        ch, dh = cs.hist.get(key, zero_ops()), ds.hist.get(key, zero_ops())
+        tot = lambda h, fld, paths=PATHS: sum(h[o][p][fld] for o in OPS for p in paths)
+        ccalls = tot(ch, "calls") - tot(ch, "calls", ("dsa_failed",))
+        dcalls = tot(dh, "calls") - tot(dh, "calls", ("dsa_failed",))
+        buckets.append({
+            "min": key[0], "max": key[1], "calls": ccalls,
+            "cpu_only_us": ratio(tot(ch, "ns"), ccalls * 1e3),
+            "with_dsa_us": ratio(tot(dh, "ns"), dcalls * 1e3),
+            "on_dsa": ratio(tot(dh, "calls", ("dsa",)), dcalls),
+            "cpu_only_s": tot(ch, "ns") / 1e9,
+        })
+    size_notes(buckets)
+    buckets.sort(key=lambda b: -b["cpu_only_s"])
+    res["buckets"] = sorted(buckets[:top], key=lambda b: b["min"])
+    return res, m
+
+
+def size_notes(buckets):
+    """Point out the sizes where DSA was slower than the CPU."""
+    compared = [b for b in buckets
+                if b["on_dsa"] and b["cpu_only_us"] and b["with_dsa_us"]]
+    # within 5% counts as neither
+    lost = [b for b in compared if b["with_dsa_us"] > b["cpu_only_us"] * 1.05]
+    won = [b for b in compared if b["cpu_only_us"] > b["with_dsa_us"] * 1.05]
+    if not lost:
+        return
+    if not won:
+        note("mem* calls that went to DSA were slower than on the CPU at every "
+             "size this program used")
+    elif max(b["min"] for b in lost) < min(b["min"] for b in won):
+        cross = min(b["min"] for b in won)
+        note(f"DSA was slower than the CPU for calls under {fmt_size(cross)} and "
+             f"faster from there up: DTO_MIN_BYTES={cross} (with "
+             "DTO_AUTO_ADJUST_KNOBS=0) keeps the smaller ones on the CPU")
+    else:
+        note("DSA was slower than the CPU for some sizes and faster for others; "
+             "see the per-size table")
+
+
+# ---- perf ------------------------------------------------------------------
+
+CATEGORIES = ("memcpy/memmove", "memset", "memcmp", "CRC", "kernel copy/clear",
+              "DTO library", "libc, unknown function", "other user code",
+              "other kernel code")
+MEM_CATEGORIES = CATEGORIES[:3]
+USER_PATTERNS = (
+    ("memcmp", re.compile(r"mem(cmp|cmpeq)|bcmp")),
+    ("memset", re.compile(r"memset|bzero")),
+    ("memcpy/memmove", re.compile(r"mem(cpy|move|pcpy)")),
+    ("CRC", re.compile(r"crc", re.I)),
+)
+KERNEL_COPY = re.compile(r"copy_(user|to_user|from_user|page|mc)|rep_movs|rep_stos|"
+                         r"clear_(page|user)|^_*mem(cpy|move|set)")
+
+
+def is_libc(path):
+    """libc.so.6, or libc-2.28.so on older glibc"""
+    return re.match(r"libc(\.so|-[0-9.]+\.so)", os.path.basename(path)) is not None
+
+
+def read_elf(path):
+    """Build ID, PT_LOAD segments and dynamic function addresses of a 64-bit
+    little-endian ELF file, or None."""
+    try:
+        with open(path, "rb") as f:
+            data = f.read()
+    except OSError:
+        return None
+    if data[:6] != b"\x7fELF\x02\x01":
+        return None
+    phoff, shoff = struct.unpack_from("<QQ", data, 0x20)
+    phentsize, phnum, shentsize, shnum = struct.unpack_from("<HHHH", data, 0x36)
+    loads, build_id, funcs = [], None, []
+    for i in range(phnum):
+        p_type, _, p_off, p_vaddr, _, p_filesz, _, _ = struct.unpack_from(
+            "<IIQQQQQQ", data, phoff + i * phentsize)
+        if p_type == 1:
+            loads.append((p_off, p_vaddr, p_filesz))
+        elif p_type == 4:
+            off, end = p_off, p_off + p_filesz
+            while off + 12 <= end:
+                namesz, descsz, ntype = struct.unpack_from("<III", data, off)
+                desc = off + 12 + ((namesz + 3) & ~3)
+                if ntype == 3 and data[off + 12:off + 15] == b"GNU":
+                    build_id = data[desc:desc + descsz].hex()
+                off = desc + ((descsz + 3) & ~3)
+    for i in range(shnum):
+        _, sh_type, _, _, sh_off, sh_size, _, _, _, _ = struct.unpack_from(
+            "<IIQQQQIIQQ", data, shoff + i * shentsize)
+        if sh_type == 11:       # SHT_DYNSYM
+            for j in range(sh_size // 24):
+                _, info, _, _, value, _ = struct.unpack_from("<IBBHQQ", data, sh_off + j * 24)
+                if (info & 0xf) in (2, 10) and value:   # STT_FUNC, STT_GNU_IFUNC
+                    funcs.append(value)
+    return {"build_id": build_id, "loads": loads, "funcs": funcs}
+
+
+class LibcLocator:
+    """Names the mem* functions in a profile of a libc without a symbol table.
+
+    Their implementations are internal (__memmove_avx512_unaligned_erms and
+    so on), so perf shows only an address. This process's libc says, through
+    dlsym, which implementation of memcpy, memset, ... this CPU uses; a sample
+    from that address up to the next function known is counted as that
+    function. Only used for a profile of the same libc file (same build ID),
+    taken on this machine: another CPU may pick other implementations."""
+
+    COPY = ("memcpy", "memmove", "mempcpy", "__memcpy_chk", "__memmove_chk",
+            "__mempcpy_chk", "wmemcpy", "wmemmove")
+    FILL = ("memset", "__memset_chk", "bzero", "wmemset", "__wmemset_chk")
+    CMP = ("memcmp", "bcmp", "__memcmpeq", "wmemcmp")
+    # other string functions, only as ends of the ones above
+    OTHER = ("strlen", "strnlen", "strchr", "strchrnul", "strrchr", "strcmp", "strncmp",
+             "strcpy", "stpcpy", "strncpy", "stpncpy", "strcat", "strncat", "memchr",
+             "memrchr", "rawmemchr", "strcasecmp", "strncasecmp", "strcspn", "strspn",
+             "strpbrk", "strstr", "wcslen", "wcsnlen", "wcschr", "wcsrchr", "wcscmp",
+             "wcsncmp", "wcscpy", "wmemchr")
+    MAX_FUNC = 64 << 10
+
+    def __init__(self):
+        self.ok = False
+        path = base = None
+        try:
+            libc = ctypes.CDLL("libc.so.6")
+            with open("/proc/self/maps") as f:
+                for line in f:
+                    fld = line.split()
+                    if len(fld) > 5 and is_libc(fld[5]) and int(fld[2], 16) == 0:
+                        path, base = fld[5], int(fld[0].split("-")[0], 16)
+                        break
+        except (OSError, ValueError):
+            return
+        elf = read_elf(path) if path else None
+        if not elf or not elf["build_id"]:
+            return
+        self.build_id, self.loads = elf["build_id"], elf["loads"]
+        starts = {a: None for a in elf["funcs"]}
+        for cat, names in (("memcpy/memmove", self.COPY), ("memset", self.FILL),
+                           ("memcmp", self.CMP), (None, self.OTHER)):
+            for name in names:
+                try:
+                    addr = ctypes.cast(getattr(libc, name), ctypes.c_void_p).value
+                except AttributeError:
+                    continue
+                starts[addr - base] = cat
+        self.addrs = sorted(starts)
+        self.cats = [starts[a] for a in self.addrs]
+        self.ok = True
+
+    def category(self, offset):
+        """Category of a sample at this file offset in libc, or None."""
+        for p_off, p_vaddr, size in self.loads:
+            if p_off <= offset < p_off + size:
+                addr = offset - p_off + p_vaddr
+                break
+        else:
+            return None
+        i = bisect.bisect_right(self.addrs, addr) - 1
+        if i < 0 or addr - self.addrs[i] > self.MAX_FUNC:
+            return None
+        return self.cats[i]
+
+
+def perf_lines(perf, args):
+    try:
+        p = subprocess.run([perf] + args, capture_output=True, text=True,
+                           errors="replace")
+    except OSError as e:
+        note(f"cannot run {perf}: {e}")
+        return None
+    if p.returncode:
+        note(f"{perf} {args[0]} failed: {p.stderr.strip()[:200]}")
+        return None
+    return p.stdout.splitlines()
+
+
+def profile(perf, data, locator):
+    """Share of samples in each category for a perf.data file, or None."""
+    libc_id = None
+    if locator.ok:
+        for line in perf_lines(perf, ["buildid-list", "-f", "-i", data]) or []:
+            fld = line.split(None, 1)
+            if len(fld) == 2 and is_libc(fld[1]):
+                libc_id = fld[0]
+    use_locator = locator.ok and libc_id == locator.build_id
+
+    lines = perf_lines(perf, ["report", "-f", "-i", data, "--stdio", "-q",
+                              "--no-children", "--sort", "dso,sym",
+                              "-F", "period,dso,sym", "-t", "|"])
+    if lines is None:
+        return None
+    shares, unknown_libc = classify(lines, locator if use_locator else None)
+    if unknown_libc > 0.05 and not use_locator:
+        note("the profiles' libc has no symbol table, so libc functions show as "
+             "addresses; run the report on the machine that ran the program (to "
+             "place them by address) or install libc debug symbols (libc6-dbg, "
+             "glibc-debuginfo)")
+    return shares
+
+
+def classify(lines, locator=None):
+    """Shares of the samples in each category, from "perf report" lines of
+    period|dso|sym, and the share in libc that could not be placed."""
+    totals = dict.fromkeys(CATEGORIES, 0)
+    unknown_libc = 0
+    for line in lines:
+        fld = line.split("|", 2)
+        if len(fld) != 3:
+            continue
+        try:
+            period = int(fld[0])
+        except ValueError:
+            continue
+        dso, sym = fld[1].strip(), fld[2].strip()
+        m = re.match(r"\[(.)\]\s*(.*)", sym)
+        kernel, name = (m.group(1) == "k", m.group(2)) if m else (False, sym)
+        if kernel:
+            cat = "kernel copy/clear" if KERNEL_COPY.search(name) else "other kernel code"
+        elif dso.startswith("libdto"):
+            cat = "DTO library"
+        else:
+            cat = next((c for c, pat in USER_PATTERNS if pat.search(name)), None)
+            if cat is None and is_libc(dso) and name.startswith("0x"):
+                cat = locator.category(int(name, 16)) if locator else None
+                if cat is None:
+                    cat = "libc, unknown function"
+                    unknown_libc += period
+            cat = cat or "other user code"
+        totals[cat] += period
+    total = sum(totals.values())
+    if not total:
+        return None, 0
+    return {k: v / total for k, v in totals.items()}, unknown_libc / total
+
+
+def profiles(d, header, runs, perf):
+    perf = perf or header.get("perf") or "perf"
+    locator = LibcLocator()
+    out = {}
+    for r in runs:
+        if r.get("profiled") and r.get("perf_data") and r.get("exit_status") == 0:
+            shares = profile(perf, os.path.join(d, r["perf_data"]), locator)
+            if shares:
+                out[r["mode"]] = shares
+    return out
+
+
+def profile_notes(res, prof):
+    cpu = prof.get("cpu")
+    if not cpu:
+        return
+    seen = res.get("amdahl", {}).get("f")
+    mem = sum(cpu[c] for c in MEM_CATEGORIES)
+    if seen is not None and mem - seen > 0.05:
+        note(f"perf puts {mem:.0%} of the CPU-only run in memcpy/memmove/memset/memcmp, "
+             f"DTO saw {seen:.0%}: the rest is calls DTO doesn't intercept, such as "
+             "calls inside libc or other libraries, __memcpy_chk and friends from "
+             "_FORTIFY_SOURCE, or a statically linked libc")
+    if cpu["CRC"] > 0.02:
+        note(f"{cpu['CRC']:.0%} of the CPU-only run is CRC code: DSA can compute "
+             "CRC32C, but DTO doesn't intercept CRC functions, so offloading it needs "
+             "an explicit call")
+    if cpu["kernel copy/clear"] > 0.02:
+        note(f"{cpu['kernel copy/clear']:.0%} of the CPU-only run is the kernel "
+             "copying or clearing memory (read/write, page faults), which DTO doesn't "
+             "offload")
+
+
+# ---- output ----------------------------------------------------------------
+
+def fmt_bytes(n):
+    for unit in ("B", "KiB", "MiB", "GiB", "TiB"):
+        if abs(n) < 1024 or unit == "TiB":
+            return f"{n:.0f} {unit}" if unit == "B" else f"{n:.1f} {unit}"
+        n /= 1024
+
+
+def fmt_size(n):
+    if n is None:
+        return ""
+    for unit in ("", "K", "M", "G"):
+        if n < 1024 or unit == "G":
+            return f"{n:g}{unit}"
+        n /= 1024
+
+
+def fmt_count(n):
+    for unit, div in (("G", 1e9), ("M", 1e6), ("k", 1e3)):
+        if n >= div:
+            return f"{n / div:.2f}{unit}"
+    return f"{n:.0f}"
+
+
+def fmt_s(s):
+    if s is None:
+        return "-"
+    if abs(s) >= 1:
+        return f"{s:.2f} s"
+    if abs(s) >= 1e-3:
+        return f"{s * 1e3:.2f} ms"
+    return f"{s * 1e6:.1f} us"
+
+
+def fmt_range(r):
+    lo, hi = fmt_s(r[0]).split(), fmt_s(r[1]).split()
+    if lo[1] == hi[1]:
+        return f"{lo[0]}-{hi[0]} {hi[1]}"
+    return f"{fmt_s(r[0])}-{fmt_s(r[1])}"
+
+
+def fmt_x(x):
+    return "-" if x is None else f"{x:.2f}x"
+
+
+def fmt_pct(x):
+    return "-" if x is None else f"{x:.1%}"
+
+
+def print_report(d, header, res, prof):
+    print(f"DTO stats run: {d}")
+    print(f"  command: {' '.join(res['command'])}")
+    print(f"  library: {res['library']}   host: {res['host']}")
+    env = " ".join(f"{k}={v}" for k, v in sorted(res["env"].items())
+                   if k not in ("DTO_COLLECT_STATS", "DTO_STATS_FORMAT"))
+    if env:
+        print(f"  settings: {env}")
+    modes = res["modes"]
+    if modes.get("cpu") and modes.get("dsa"):
+        print(f"  timed runs: {modes['cpu']['runs']} CPU-only, {modes['dsa']['runs']} "
+              f"with DSA (medians)")
+    cfg = res.get("dsa_config")
+    if cfg:
+        print(f"  DTO in the DSA run: {cfg.get('wqs', '?')} WQs, wait "
+              f"{cfg.get('wait_method', '?')}, min_bytes {cfg.get('min_bytes', '?')}, "
+              f"cpu_size_fraction {cfg.get('cpu_size_fraction', '?')}"
+              f"{' (auto-adjusted)' if cfg.get('auto_adjust') else ''}")
+    print()
+
+    if "amdahl" in res:
+        c, dm = modes["cpu"], modes["dsa"]
+        out = res["outside_mem_s"]
+        print(f"  {'':28}{'CPU only':>15}{'with DSA':>15}{'speedup':>10}")
+        print(f"  {'wall time':28}{fmt_s(c['wall_s']):>15}{fmt_s(dm['wall_s']):>15}"
+              f"{fmt_x(res['measured']['wall']):>10}")
+        print(f"  {'CPU time (user+sys)':28}{fmt_s(c['cpu_s']):>15}{fmt_s(dm['cpu_s']):>15}"
+              f"{fmt_x(res['measured']['cpu']):>10}")
+        print(f"  {'  in mem* calls':28}{fmt_s(c['mem_s']):>15}{fmt_s(dm['mem_s']):>15}"
+              f"{fmt_x(res['amdahl']['s']):>10}")
+        print(f"  {'  outside mem* calls':28}{fmt_s(out['cpu']):>15}{fmt_s(out['dsa']):>15}"
+              f"{fmt_x(ratio(out['cpu'], out['dsa']) if out['dsa'] > 0 else None):>10}")
+        if c["runs"] > 1 or dm["runs"] > 1:
+            print(f"  {'wall time, min-max':28}{fmt_range(c['wall_range']):>15}"
+                  f"{fmt_range(dm['wall_range']):>15}")
+        a = res["amdahl"]
+        print()
+        print("  Amdahl's law")
+        print(f"    f, share of CPU time in mem* calls (CPU-only)  {fmt_pct(a['f']):>8}")
+        print(f"    s, mem* speedup with DSA (time per byte)       {fmt_x(a['s']):>8}")
+        print(f"    predicted speedup, 1 / ((1 - f) + f / s)       {fmt_x(a['predicted']):>8}")
+        print(f"    measured speedup, wall time                    "
+              f"{fmt_x(res['measured']['wall']):>8}")
+        print(f"    limit if mem* took no time, 1 / (1 - f)        {fmt_x(a['limit']):>8}")
+
+    if res.get("ops"):
+        print()
+        print(f"  {'per operation':16}{'calls':>9}{'bytes':>12}{'avg size':>11}"
+              f"{'CPU only':>11}{'with DSA':>11}{'speedup':>9}{'on DSA':>8}{'f':>8}")
+        for op, o in res["ops"].items():
+            print(f"  {op:16}{fmt_count(o['calls']):>9}{fmt_bytes(o['bytes']):>12}"
+                  f"{fmt_bytes(o['bytes'] / o['calls']) if o['calls'] else '-':>11}"
+                  f"{fmt_s(o['cpu_only_s']):>11}{fmt_s(o['with_dsa_s']):>11}"
+                  f"{fmt_x(o['speedup']):>9}{fmt_pct(o['on_dsa']):>8}{fmt_pct(o['share']):>8}")
+        print("    (on DSA: share of the bytes in calls sent to DSA, of which the CPU "
+              "does cpu_size_fraction;")
+        print("     f: share of the CPU-only run's CPU time)")
+        fails = res.get("dsa_failures", {})
+        if any(fails.values()):
+            print(f"    calls DSA left to the CPU, per run: page fault "
+                  f"{fmt_count(fails['page_fault'])}, retry {fmt_count(fails['retry'])}, "
+                  f"other {fmt_count(fails['other'])}")
+
+    if res.get("buckets"):
+        print()
+        print(f"  {'by size (most time)':22}{'calls':>9}{'CPU only':>11}{'with DSA':>11}"
+              f"{'speedup':>9}{'on DSA':>8}")
+        for b in res["buckets"]:
+            # buckets are [min, max]; [min, max + 1) reads better
+            rng = f"{fmt_size(b['min'])}-{fmt_size(b['max'] + 1)}" if b["max"] else \
+                f">={fmt_size(b['min'])}"
+            us = lambda v: "-" if v is None else f"{v:.2f} us"
+            print(f"  {rng:22}{fmt_count(b['calls']):>9}{us(b['cpu_only_us']):>11}"
+                  f"{us(b['with_dsa_us']):>11}"
+                  f"{fmt_x(ratio(b['cpu_only_us'], b['with_dsa_us'])):>9}"
+                  f"{fmt_pct(b['on_dsa']):>8}")
+        print("    (average time per call; on DSA: share of calls sent to DSA)")
+
+    if prof:
+        print()
+        print(f"  {'profile (perf, % of CPU time)':34}"
+              + "".join(f"{('CPU only' if m == 'cpu' else 'with DSA'):>11}"
+                        for m in MODES if m in prof))
+        for cat in CATEGORIES:
+            if not any(prof[m][cat] >= 0.0005 for m in prof):
+                continue
+            print(f"    {cat:32}" + "".join(f"{prof[m][cat]:>11.1%}"
+                                             for m in MODES if m in prof))
+
+    if notes:
+        print()
+        for n in notes:
+            print(f"  note: {n}")
+
+
+def main():
+    ap = argparse.ArgumentParser(
+        description="Summarize a dtoctl --stats-run: DSA speedup and Amdahl's law.")
+    ap.add_argument("dir", help="directory written by dtoctl --stats-run")
+    ap.add_argument("--json", action="store_true", help="print the results as JSON")
+    ap.add_argument("--perf", help="perf binary for reading the profiles "
+                    "(default: the one dtoctl ran)")
+    ap.add_argument("--top", type=int, default=8,
+                    help="size ranges to show, those with the most CPU-only time "
+                    "(default: 8)")
+    args = ap.parse_args()
+
+    header, runs = load_dir(args.dir)
+    res, _ = analyze(header, runs, args.top)
+    prof = profiles(args.dir, header, runs, args.perf)
+    profile_notes(res, prof)
+    if args.json:
+        res["profile"] = prof
+        res["notes"] = notes
+        json.dump(res, sys.stdout, indent=1)
+        print()
+    else:
+        print_report(args.dir, header, res, prof)
+
+
+if __name__ == "__main__":
+    main()
diff --git a/dto.c b/dto.c
index dd76ef4..8c7f4dd 100644
--- a/dto.c
+++ b/dto.c
@@ -4,11 +4,13 @@
  * SPDX-License-Identifier: MIT
  ******************************************************************************/
 
+#include <stddef.h>
 #include <stdint.h>
 #include <stdio.h>
 #include <stdlib.h>
 #include <fcntl.h>
 #include <errno.h>
+#include <poll.h>
 #include <string.h>
 #include <stdarg.h>
 #include <unistd.h>
@@ -187,9 +189,13 @@ static const char * const memop_names[] = {
 	[MEMCMP] = "cmp"
 };
 
-// memory stats
-#define HIST_BUCKET_SIZE 4096
-#define HIST_NO_BUCKETS 512
+// memory stats: a histogram by call size, with bucket 0 for calls under
+// 64 bytes, then four buckets per power of two up to 4 GiB, and a last
+// bucket for anything larger
+#define HIST_MIN_SHIFT 6
+#define HIST_MAX_SHIFT 32
+#define HIST_SUBS 4
+#define HIST_NO_BUCKETS (2 + (HIST_MAX_SHIFT - HIST_MIN_SHIFT) * HIST_SUBS)
 enum stat_group {
 	STDC_CALL = 0x0,
 	DSA_CALL_SUCCESS,
@@ -268,10 +274,58 @@ static struct timespec dto_start_time;
 		}								\
 	} while (0)								\
 
-static atomic_int op_counter[HIST_NO_BUCKETS][MAX_STAT_GROUP][MAX_MEMOP];
-static atomic_ullong bytes_counter[HIST_NO_BUCKETS][MAX_STAT_GROUP];
-static atomic_ullong lat_counter[HIST_NO_BUCKETS][MAX_STAT_GROUP][MAX_MEMOP];
-static atomic_int fail_counter[HIST_NO_BUCKETS][MAX_FAILURES];
+/* the groups that have per-operation counters */
+#define STAT_GROUPS DSA_FAIL_CODES
+
+struct dto_op_stats {
+	uint64_t calls;
+	uint64_t bytes;
+	uint64_t ns;
+};
+
+struct dto_bucket_stats {
+	struct dto_op_stats op[STAT_GROUPS][MAX_MEMOP];
+	uint64_t fails[MAX_FAILURES];
+};
+
+struct dto_stats {
+	struct dto_op_stats total[STAT_GROUPS][MAX_MEMOP];
+	uint64_t fails[MAX_FAILURES];
+	struct dto_bucket_stats hist[HIST_NO_BUCKETS];
+};
+
+/* one per counting thread (see update_stats) */
+struct dto_thread_stats {
+	struct dto_thread_stats *prev, *next;
+	struct dto_stats s;
+};
+
+static struct dto_thread_stats *stats_threads;	/* live threads' blocks */
+static struct dto_stats stats_retired;		/* exited threads, summed */
+static pthread_mutex_t stats_lock = PTHREAD_MUTEX_INITIALIZER;
+static pthread_key_t stats_key;
+static bool stats_key_created;
+static __thread struct dto_thread_stats *thr_stats;
+static __thread bool thr_stats_off;	/* don't count this thread's calls */
+
+/* output: one formatter at a time (sampler, exit, SIGINT) */
+static pthread_mutex_t stats_out_lock = PTHREAD_MUTEX_INITIALIZER;
+static struct dto_stats stats_scratch;
+enum { STATS_FORMAT_TEXT, STATS_FORMAT_JSON };
+static int stats_format;
+
+static bool stats_file_per_pid;	/* DTO_STATS_FILE has a %p */
+
+/* DTO_STATS_INTERVAL_MS */
+static unsigned long stats_interval_ms;
+static atomic_uchar sampler_started;
+static bool stats_shutdown;	/* cleanup_dto ran: no new sampler */
+static struct dto_stats sampler_prev;	/* totals at the previous record */
+static uint64_t sampler_prev_ns;
+
+static void print_stats(void);
+static void print_stats_from_signal(void);
+static void stats_fork_child(void);
 #endif
 
 /* call initialize/cleanup functions when library is loaded/unloaded */
@@ -292,29 +346,7 @@ static enum {
 
 static unsigned int log_level = LOG_LEVEL_FATAL;
 
-static void dto_stats_log(const char *fmt, ...)
-{
-        char buf[512];
-        va_list args;
-
-        if (stats_fd == -1 && LOG_LEVEL_TRACE > log_level)
-                return;
-
-        va_start(args, fmt);
-        if (stats_fd != -1) {
-                if (vsnprintf(buf, sizeof(buf), fmt, args) > 0)
-                        write(stats_fd, buf, strlen(buf));
-        } else if (log_fd == -1) {
-                vprintf(fmt, args);
-        } else {
-                vsnprintf(buf, sizeof(buf), fmt, args);
-                write(log_fd, buf, strlen(buf));
-        }
-        va_end(args);
-}
-
 #ifdef DTO_STATS_SUPPORT
-static void print_stats(void);
 static struct sigaction previous_sigint_action;
 static volatile sig_atomic_t sigint_handler_installed;
 static int previous_sigint_action_valid;
@@ -338,8 +370,7 @@ static void restore_sigint_handler(void)
 
 static void dto_handle_sigint(int signo)
 {
-        if (collect_stats)
-                print_stats();
+        print_stats_from_signal();
 
         restore_sigint_handler();
 
@@ -418,20 +449,7 @@ static void dto_log(int req_log_level, const char *fmt, ...)
 static void child (void)
 {
 #ifdef DTO_STATS_SUPPORT
-	int i, j, k;
-
-	/* Reset the counters */
-	for (i = 0; i < HIST_NO_BUCKETS; i++) {
-		for (j = 0; j < MAX_STAT_GROUP; j++) {
-			for (k = 0; k < MAX_MEMOP; k++) {
-				op_counter[i][j][k] = 0;
-				lat_counter[i][j][k] = 0;
-			}
-			bytes_counter[i][j] = 0;
-		}
-		for (j = 0; j < MAX_FAILURES; j++)
-			fail_counter[i][j] = 0;
-	}
+	stats_fork_child();
 #endif
 	dto_initializing = 0;
 	dto_initialized = 0;
@@ -672,114 +690,677 @@ static __always_inline int dsa_execute(struct dto_wq *wq,
 }
 
 #ifdef DTO_STATS_SUPPORT
+/*
+ * Statistics (DTO_COLLECT_STATS).
+ *
+ * Each thread counts into its own block, so an update is a few plain loads
+ * and stores with no shared cache line; the blocks are summed when stats are
+ * printed. A block is mmap'd on the thread's first update (only the pages it
+ * touches are ever allocated), registered in a list, and on thread exit
+ * folded into stats_retired and unmapped. The totals sit at the start of the
+ * block, so the interval sampler reads one page per thread.
+ */
+
+/* Histogram bucket of a call of n bytes. */
+static __always_inline int stats_bucket(size_t n)
+{
+	int e;
+
+	if (n < (1UL << HIST_MIN_SHIFT))
+		return 0;
+	e = 63 - __builtin_clzl(n);
+	if (e >= HIST_MAX_SHIFT)
+		return HIST_NO_BUCKETS - 1;
+	return 1 + (e - HIST_MIN_SHIFT) * HIST_SUBS +
+		(int)((n >> (e - 2)) & (HIST_SUBS - 1));
+}
+
+/* Smallest size counted in bucket b. */
+static uint64_t stats_bucket_min(int b)
+{
+	int e;
+
+	if (b == 0)
+		return 0;
+	if (b == HIST_NO_BUCKETS - 1)
+		return 1ULL << HIST_MAX_SHIFT;
+	e = HIST_MIN_SHIFT + (b - 1) / HIST_SUBS;
+	return (uint64_t)(HIST_SUBS + (b - 1) % HIST_SUBS) << (e - 2);
+}
+
+/* Relaxed load-add-store: the owning thread is the only writer, and
+ * readers summing blocks only need each counter to be read whole. */
+#define STAT_ADD(field, v) \
+	__atomic_store_n(&(field), __atomic_load_n(&(field), __ATOMIC_RELAXED) + (v), \
+			 __ATOMIC_RELAXED)
+#define STAT_GET(field) __atomic_load_n(&(field), __ATOMIC_RELAXED)
+
+static __always_inline void stats_add_op(struct dto_op_stats *s, uint64_t bytes,
+					 uint64_t ns)
+{
+	STAT_ADD(s->calls, 1);
+	STAT_ADD(s->bytes, bytes);
+	STAT_ADD(s->ns, ns);
+}
+
+/* dst += src for the totals, or for everything with hist */
+static void stats_sum(struct dto_stats *dst, struct dto_stats *src, bool hist)
+{
+	for (int g = 0; g < STAT_GROUPS; g++)
+		for (int o = 0; o < MAX_MEMOP; o++) {
+			dst->total[g][o].calls += STAT_GET(src->total[g][o].calls);
+			dst->total[g][o].bytes += STAT_GET(src->total[g][o].bytes);
+			dst->total[g][o].ns += STAT_GET(src->total[g][o].ns);
+		}
+	for (int f = 0; f < MAX_FAILURES; f++)
+		dst->fails[f] += STAT_GET(src->fails[f]);
+	if (!hist)
+		return;
+	for (int b = 0; b < HIST_NO_BUCKETS; b++) {
+		struct dto_bucket_stats *d = &dst->hist[b], *s = &src->hist[b];
+
+		for (int g = 0; g < STAT_GROUPS; g++)
+			for (int o = 0; o < MAX_MEMOP; o++) {
+				d->op[g][o].calls += STAT_GET(s->op[g][o].calls);
+				d->op[g][o].bytes += STAT_GET(s->op[g][o].bytes);
+				d->op[g][o].ns += STAT_GET(s->op[g][o].ns);
+			}
+		for (int f = 0; f < MAX_FAILURES; f++)
+			d->fails[f] += STAT_GET(s->fails[f]);
+	}
+}
+
+/* Sum of every thread's stats, live and exited, into *out (zeroed first).
+ * The caller holds stats_lock. */
+static void stats_collect(struct dto_stats *out, bool hist)
+{
+	size_t len = hist ? sizeof(*out) : offsetof(struct dto_stats, hist);
+
+	orig_memset(out, 0, len);
+	stats_sum(out, &stats_retired, hist);
+	for (struct dto_thread_stats *t = stats_threads; t; t = t->next)
+		stats_sum(out, &t->s, hist);
+}
+
+static void stats_snapshot(struct dto_stats *out, bool hist)
+{
+	pthread_mutex_lock(&stats_lock);
+	stats_collect(out, hist);
+	pthread_mutex_unlock(&stats_lock);
+}
+
+static void stats_start_sampler(void);
+
+static void stats_thread_exit(void *arg)
+{
+	struct dto_thread_stats *ts = arg;
+
+	/* memcpy and friends may still run during the rest of thread exit */
+	thr_stats_off = true;
+	thr_stats = NULL;
+	pthread_mutex_lock(&stats_lock);
+	stats_sum(&stats_retired, &ts->s, true);
+	if (ts->prev)
+		ts->prev->next = ts->next;
+	else
+		stats_threads = ts->next;
+	if (ts->next)
+		ts->next->prev = ts->prev;
+	pthread_mutex_unlock(&stats_lock);
+	munmap(ts, sizeof(*ts));
+}
+
+/* First update on this thread: give it a block. NULL leaves the thread
+ * uncounted (thr_stats_off stays set, so it is not retried every call). */
+static struct dto_thread_stats *stats_thread_register(void)
+{
+	struct dto_thread_stats *ts;
+
+	/* mmap, the registry and pthread_setspecific must not count themselves */
+	thr_stats_off = true;
+	ts = mmap(NULL, sizeof(*ts), PROT_READ | PROT_WRITE,
+		  MAP_PRIVATE | MAP_ANONYMOUS, -1, 0);
+	if (ts == MAP_FAILED)
+		return NULL;
+
+	pthread_mutex_lock(&stats_lock);
+	ts->prev = NULL;
+	ts->next = stats_threads;
+	if (stats_threads)
+		stats_threads->prev = ts;
+	stats_threads = ts;
+	pthread_mutex_unlock(&stats_lock);
+	pthread_setspecific(stats_key, ts);
+	thr_stats = ts;
+
+	if (stats_interval_ms && !stats_shutdown &&
+	    !atomic_load_explicit(&sampler_started,
+						       memory_order_acquire))
+		stats_start_sampler();
+	thr_stats_off = false;
+	return ts;
+}
+
 static void update_stats(int op, size_t n, bool overlapping, size_t bytes_completed,
 		uint64_t elapsed_ns, int group, int error_code)
 {
+	struct dto_thread_stats *ts;
+	struct dto_bucket_stats *b;
+
 	// dto_memcpymove didn't actually submit the request to DSA, so there is nothing to log. This will be captured by a second call
 	if (op == MEMMOVE && overlapping && dto_overlapping_memmove_action == OVERLAPPING_CPU && group == DSA_CALL_SUCCESS) {
 		return;
 	}
+	if (unlikely(thr_stats_off))
+		return;
+	ts = thr_stats;
+	if (unlikely(ts == NULL)) {
+		ts = stats_thread_register();
+		if (ts == NULL)
+			return;
+	}
 
-	int bucket = (n / HIST_BUCKET_SIZE);
+	b = &ts->s.hist[stats_bucket(n)];
+	stats_add_op(&ts->s.total[group][op], bytes_completed, elapsed_ns);
+	stats_add_op(&b->op[group][op], bytes_completed, elapsed_ns);
+	if (group == DSA_CALL_FAILED) {
+		STAT_ADD(ts->s.fails[error_code], 1);
+		STAT_ADD(b->fails[error_code], 1);
+	}
+}
 
-	if (bucket >= HIST_NO_BUCKETS)  /* last bucket includes remaining sizes */
-		bucket = HIST_NO_BUCKETS-1;
-	++op_counter[bucket][group][op];
-	bytes_counter[bucket][group] += bytes_completed;
-	lat_counter[bucket][group][op] += elapsed_ns;
-	if (group == DSA_CALL_FAILED)
-		++fail_counter[bucket][error_code];
+/* ---- Output ----
+ *
+ * Records are formatted into stats_out and written with as few write()s as
+ * possible: to DTO_STATS_FILE (opened O_APPEND, so records from several
+ * processes sharing it stay whole), else DTO_LOG_FILE, else stdout. */
+
+/* Holds a whole record, so a summary (up to ~1 KB for each of the 106
+ * histogram buckets in JSON) goes out in one write() and records from
+ * several processes sharing an O_APPEND file stay whole. */
+static char stats_out[1 << 18];
+static size_t stats_out_len;
+
+static void stats_flush(void)
+{
+	int fd = stats_fd != -1 ? stats_fd : log_fd != -1 ? log_fd : STDOUT_FILENO;
+	size_t off = 0;
+
+	while (off < stats_out_len) {
+		ssize_t r = write(fd, stats_out + off, stats_out_len - off);
 
+		if (r < 0 && errno == EINTR)
+			continue;
+		if (r < 0 && errno == EAGAIN) {
+			struct pollfd pfd = { .fd = fd, .events = POLLOUT };
+
+			if (poll(&pfd, 1, -1) > 0 || errno == EINTR)
+				continue;
+		}
+		if (r <= 0)
+			break;
+		off += r;
+	}
+	stats_out_len = 0;
 }
 
-static void print_stats(void)
+__attribute__((format(printf, 1, 2)))
+static void stats_printf(const char *fmt, ...)
 {
-	struct timespec dto_end_time;
+	va_list args;
+	int len;
 
-	if (likely(!collect_stats))
+	va_start(args, fmt);
+	len = vsnprintf(stats_out + stats_out_len, sizeof(stats_out) - stats_out_len,
+			fmt, args);
+	va_end(args);
+	if (len < 0)
 		return;
+	if ((size_t)len >= sizeof(stats_out) - stats_out_len) {
+		/* no room: write out what is there and format again */
+		stats_flush();
+		va_start(args, fmt);
+		len = vsnprintf(stats_out, sizeof(stats_out), fmt, args);
+		va_end(args);
+		if (len < 0)
+			return;
+		if ((size_t)len >= sizeof(stats_out))
+			len = sizeof(stats_out) - 1;
+	}
+	stats_out_len += len;
+}
+
+/* DTO_STATS_FILE, with each %p replaced by the process ID. It is opened
+ * O_APPEND, not truncated: processes started from this one (fork, exec)
+ * inherit the setting, and with %p each gets a file of its own, without it
+ * they append whole records to a shared one. */
+static void stats_open_file(const char *name)
+{
+	char path[PATH_MAX];
+	size_t j = 0;
+	struct stat st;
+
+	for (const char *c = name; *c && j < sizeof(path) - 1; c++) {
+		if (c[0] == '%' && c[1] == 'p') {
+			stats_file_per_pid = true;
+			int n = snprintf(path + j, sizeof(path) - j, "%d", getpid());
+
+			if (n < 0 || (size_t)n >= sizeof(path) - j)
+				break;
+			j += n;
+			c++;
+		} else {
+			path[j++] = *c;
+		}
+	}
+	path[j] = '\0';
+
+	if (lstat(path, &st) == 0 && (st.st_mode & S_IFMT) != S_IFREG) {
+		LOG_ERROR("DTO_STATS_FILE '%s' must reference a regular file\n", path);
+		return;
+	}
+	stats_fd = open(path, O_WRONLY | O_CREAT | O_APPEND | O_CLOEXEC, 0600);
+	if (stats_fd == -1)
+		LOG_ERROR("Failed to open DTO_STATS_FILE '%s': %s\n", path, strerror(errno));
+}
+
+static const char * const memop_json[] = {
+	[MEMSET] = "memset",
+	[MEMCOPY] = "memcpy",
+	[MEMMOVE] = "memmove",
+	[MEMCMP] = "memcmp",
+};
+
+static const char * const group_json[] = {
+	[STDC_CALL] = "cpu",
+	[DSA_CALL_SUCCESS] = "dsa",
+	[DSA_CALL_FAILED] = "dsa_failed",
+};
+
+static const char * const failure_json[] = {
+	[RETRY] = "retry",
+	[PAGE_FAULT] = "page_fault",
+	[FAIL_OTHERS] = "other",
+};
+
+static uint64_t stats_now_ns(void)
+{
+	struct timespec now;
+
+	clock_gettime(CLOCK_BOOTTIME, &now);
+	return TS_NS(dto_start_time, now);
+}
+
+/* {"memcpy":{"dsa":{"calls":..,"bytes":..,"ns":..},..},..} without zeroes */
+static void stats_json_ops(struct dto_op_stats op[STAT_GROUPS][MAX_MEMOP])
+{
+	bool first_op = true;
+
+	stats_printf("{");
+	for (int o = 0; o < MAX_MEMOP; o++) {
+		bool first_group = true;
+
+		for (int g = 0; g < STAT_GROUPS; g++) {
+			if (!op[g][o].calls)
+				continue;
+			if (first_group)
+				stats_printf("%s\"%s\":{", first_op ? "" : ",", memop_json[o]);
+			stats_printf("%s\"%s\":{\"calls\":%lu,\"bytes\":%lu,\"ns\":%lu}",
+				     first_group ? "" : ",", group_json[g], op[g][o].calls,
+				     op[g][o].bytes, op[g][o].ns);
+			first_group = false;
+			first_op = false;
+		}
+		if (!first_group)
+			stats_printf("}");
+	}
+	stats_printf("}");
+}
+
+static void stats_json_fails(uint64_t fails[MAX_FAILURES])
+{
+	stats_printf("{");
+	for (int f = 1; f < MAX_FAILURES; f++)
+		stats_printf("%s\"%s\":%lu", f == 1 ? "" : ",", failure_json[f], fails[f]);
+	stats_printf("}");
+}
 
-	clock_gettime(CLOCK_BOOTTIME, &dto_end_time);
+/* Text: "memcpy dsa 120 calls 7864320 B 421.3 us, cpu ..." */
+static void stats_text_ops(struct dto_op_stats op[STAT_GROUPS][MAX_MEMOP])
+{
+	bool any = false;
 
-	dto_stats_log("DTO Run Time: %ld ms\n", TS_NS(dto_start_time, dto_end_time)/1000000);
+	for (int o = 0; o < MAX_MEMOP; o++)
+		for (int g = 0; g < STAT_GROUPS; g++) {
+			if (!op[g][o].calls)
+				continue;
+			stats_printf("%s%s %s %lu calls %lu B %.1f us", any ? ", " : "",
+				     memop_json[o], group_json[g], op[g][o].calls,
+				     op[g][o].bytes, op[g][o].ns / 1000.0);
+			any = true;
+		}
+	if (!any)
+		stats_printf("no calls");
+}
+
+/* One interval record: the change in the totals since the previous one.
+ * The caller holds stats_out_lock. */
+static void stats_emit_interval(void)
+{
+	struct dto_stats *now = &stats_scratch;	/* totals part only */
+	struct dto_op_stats delta[STAT_GROUPS][MAX_MEMOP];
+	uint64_t fails[MAX_FAILURES];
+	uint64_t t = stats_now_ns();
+
+	stats_snapshot(now, false);
+	for (int g = 0; g < STAT_GROUPS; g++)
+		for (int o = 0; o < MAX_MEMOP; o++) {
+			delta[g][o].calls = now->total[g][o].calls - sampler_prev.total[g][o].calls;
+			delta[g][o].bytes = now->total[g][o].bytes - sampler_prev.total[g][o].bytes;
+			delta[g][o].ns = now->total[g][o].ns - sampler_prev.total[g][o].ns;
+		}
+	for (int f = 0; f < MAX_FAILURES; f++)
+		fails[f] = now->fails[f] - sampler_prev.fails[f];
+	orig_memcpy(&sampler_prev, now, offsetof(struct dto_stats, hist));
+
+	if (stats_format == STATS_FORMAT_JSON) {
+		stats_printf("{\"type\":\"interval\",\"pid\":%d,\"t_ms\":%.3f,\"dt_ms\":%.3f,\"ops\":",
+			     getpid(), t / 1e6, (t - sampler_prev_ns) / 1e6);
+		stats_json_ops(delta);
+		stats_printf(",\"failures\":");
+		stats_json_fails(fails);
+		stats_printf("}\n");
+	} else {
+		stats_printf("DTO stats t=%.3f ms (+%.3f ms): ", t / 1e6, (t - sampler_prev_ns) / 1e6);
+		stats_text_ops(delta);
+		if (fails[RETRY] || fails[PAGE_FAULT] || fails[FAIL_OTHERS])
+			stats_printf("; failures retry %lu page_fault %lu other %lu",
+				     fails[RETRY], fails[PAGE_FAULT], fails[FAIL_OTHERS]);
+		stats_printf("\n");
+	}
+	sampler_prev_ns = t;
+	stats_flush();
+}
+
+static bool stats_bucket_empty(struct dto_bucket_stats *b)
+{
+	for (int g = 0; g < STAT_GROUPS; g++)
+		for (int o = 0; o < MAX_MEMOP; o++)
+			if (b->op[g][o].calls)
+				return false;
+	return true;
+}
+
+static void stats_print_json(struct dto_stats *s, uint64_t run_ns)
+{
+	char prog[64];
+	bool first = true;
+	size_t j = 0;
+
+	/* the program name as a JSON string: drop quotes, backslashes and
+	 * control characters rather than escape them */
+	for (const char *c = __progname; *c && j < sizeof(prog) - 1; c++)
+		if (*c != '"' && *c != '\\' && (unsigned char)*c >= 0x20)
+			prog[j++] = *c;
+	prog[j] = '\0';
+	stats_printf("{\"type\":\"summary\",\"pid\":%d,\"program\":\"%s\",\"run_ms\":%.3f,",
+		     getpid(), prog, run_ns / 1e6);
+	stats_printf("\"config\":{\"stdc_only\":%d,\"wait_method\":\"%s\",\"min_bytes\":%lu,"
+		     "\"cpu_size_fraction\":%.2f,\"auto_adjust\":%d,\"numa_aware\":%d,\"wqs\":%d},",
+		     use_std_lib_calls, wait_names[wait_method], dsa_min_size,
+		     cpu_size_fraction / 100.0, auto_adjust_knobs, (int)is_numa_aware, num_wqs);
+	stats_printf("\"totals\":");
+	stats_json_ops(s->total);
+	stats_printf(",\"failures\":");
+	stats_json_fails(s->fails);
+	stats_printf(",\"histogram\":[");
+	for (int b = 0; b < HIST_NO_BUCKETS; b++) {
+		if (stats_bucket_empty(&s->hist[b]))
+			continue;
+		stats_printf("%s{\"min\":%lu,", first ? "" : ",", stats_bucket_min(b));
+		if (b < HIST_NO_BUCKETS - 1)
+			stats_printf("\"max\":%lu,", stats_bucket_min(b + 1) - 1);
+		stats_printf("\"ops\":");
+		stats_json_ops(s->hist[b].op);
+		stats_printf(",\"failures\":");
+		stats_json_fails(s->hist[b].fails);
+		stats_printf("}");
+		first = false;
+	}
+	stats_printf("]}\n");
+}
+
+static void stats_print_text(struct dto_stats *s, uint64_t run_ns)
+{
+	stats_printf("DTO Run Time: %lu ms\n", run_ns / 1000000);
 
-	// display stats
 	for (int t = 0; t < 2; ++t) {
 		if (t == 0)
-			dto_stats_log("\n******** Number of Memory Operations ********\n");
+			stats_printf("\n******** Number of Memory Operations ********\n");
 		else
-			dto_stats_log("\n******** Average Memory Operation Latency (us)  ********\n");
+			stats_printf("\n******** Average Memory Operation Latency (us)  ********\n");
 
-		dto_stats_log("%17s    ", "");
+		stats_printf("%23s    ", "");
 		for (int g = 0; g < MAX_STAT_GROUP; ++g) {
 			if (t == 0) {
 				if (g == DSA_FAIL_CODES)
-					dto_stats_log("<***** %-13s *****> ", stat_group_names[g]);
+					stats_printf("<***** %-13s *****> ", stat_group_names[g]);
 				else
-					dto_stats_log("<*************** %-13s ***************> ", stat_group_names[g]);
-			} else {
-				if (g != DSA_FAIL_CODES)
-					dto_stats_log("<******** %-13s ********> ", stat_group_names[g]);
-
+					stats_printf("<*************** %-13s ***************> ", stat_group_names[g]);
+			} else if (g != DSA_FAIL_CODES) {
+				stats_printf("<******** %-13s ********> ", stat_group_names[g]);
 			}
 		}
-		dto_stats_log("\n");
-
-		dto_stats_log("%-17s -- ", "Byte Range");
-		for (int g = 0; g < MAX_STAT_GROUP - 1; ++g) {
+		stats_printf("\n%-23s -- ", "Byte Range");
+		for (int g = 0; g < STAT_GROUPS; ++g) {
 			for (int o = 0; o < MAX_MEMOP; ++o)
-				dto_stats_log("%-8s ", memop_names[o]);
+				stats_printf("%-8s ", memop_names[o]);
 			if (t == 0)
-				dto_stats_log("%-12s ", "bytes");
+				stats_printf("%-12s ", "bytes");
 		}
 		if (t == 0)
 			for (int o = 1; o < MAX_FAILURES; ++o)
-				dto_stats_log("%-6s ", failure_names[o]);
-		dto_stats_log("\n");
+				stats_printf("%-6s ", failure_names[o]);
+		stats_printf("\n");
 
 		for (int b = 0; b < HIST_NO_BUCKETS; ++b) {
-			bool empty = true;
+			struct dto_bucket_stats *h = &s->hist[b];
 
-			for (int g = 0; g < MAX_STAT_GROUP; ++g) {
-				for (int o = 0; o < MAX_MEMOP; ++o) {
-					if (op_counter[b][g][o] != 0) {
-						empty = false;
-						break;
-					}
-				}
-				if (!empty)
-					break;
-			}
-			if (empty)
+			if (stats_bucket_empty(h))
 				continue;
-
-			if (b < (HIST_NO_BUCKETS-1))
-				dto_stats_log("% 8d-%-8d -- ", b*4096, ((b+1)*4096)-1);
+			if (b < HIST_NO_BUCKETS - 1)
+				stats_printf("%11lu-%-11lu -- ", stats_bucket_min(b),
+					     stats_bucket_min(b + 1) - 1);
 			else
-				dto_stats_log("   >=%-12d -- ", b*4096);
+				stats_printf("   >=%-18lu -- ", stats_bucket_min(b));
+
+			for (int g = 0; g < STAT_GROUPS; ++g) {
+				uint64_t bytes = 0;
 
-			for (int g = 0; g < MAX_STAT_GROUP - 1; ++g) {
 				for (int o = 0; o < MAX_MEMOP; ++o) {
-					if (t == 0) {
-						dto_stats_log("%-8d ", op_counter[b][g][o]);
-						continue;
-					}
-					if (op_counter[b][g][o] != 0) {
-						double avg_us = ((double) lat_counter[b][g][o])/(((double) op_counter[b][g][o]) * 1000.0);
+					struct dto_op_stats *op = &h->op[g][o];
 
-						dto_stats_log("%-8.2f ", avg_us);
-					} else {
-						dto_stats_log("%-8d ", 0);
-					}
+					bytes += op->bytes;
+					if (t == 0)
+						stats_printf("%-8lu ", op->calls);
+					else if (op->calls)
+						stats_printf("%-8.2f ", op->ns / (op->calls * 1000.0));
+					else
+						stats_printf("%-8d ", 0);
 				}
 				if (t == 0)
-					dto_stats_log("%-12lld ", bytes_counter[b][g]);
+					stats_printf("%-12lu ", bytes);
 			}
 			if (t == 0)
 				for (int o = 1; o < MAX_FAILURES; ++o)
-					dto_stats_log("%-6d ", fail_counter[b][o]);
-			dto_stats_log("\n");
+					stats_printf("%-6lu ", h->fails[o]);
+			stats_printf("\n");
+		}
+	}
+}
+
+/* The end-of-run summary of *s. The caller holds stats_out_lock. */
+static void stats_print_summary(struct dto_stats *s)
+{
+	if (stats_format == STATS_FORMAT_JSON)
+		stats_print_json(s, stats_now_ns());
+	else
+		stats_print_text(s, stats_now_ns());
+	stats_flush();
+}
+
+static void print_stats(void)
+{
+	bool off = thr_stats_off;
+
+	if (likely(!collect_stats))
+		return;
+	thr_stats_off = true;
+	pthread_mutex_lock(&stats_out_lock);
+	stats_snapshot(&stats_scratch, true);
+	stats_print_summary(&stats_scratch);
+	pthread_mutex_unlock(&stats_out_lock);
+	thr_stats_off = off;
+}
+
+/* From the SIGINT handler: the interrupted thread may hold a stats lock, so
+ * print only if both can be taken right away. */
+static void print_stats_from_signal(void)
+{
+	bool off = thr_stats_off;
+
+	if (!collect_stats)
+		return;
+	thr_stats_off = true;
+	if (pthread_mutex_trylock(&stats_out_lock) == 0) {
+		if (pthread_mutex_trylock(&stats_lock) == 0) {
+			stats_collect(&stats_scratch, true);
+			pthread_mutex_unlock(&stats_lock);
+			stats_print_summary(&stats_scratch);
 		}
+		pthread_mutex_unlock(&stats_out_lock);
+	}
+	thr_stats_off = off;
+}
+
+/* ---- Interval sampler (DTO_STATS_INTERVAL_MS) ---- */
+
+static pthread_t sampler_thread;
+static pthread_mutex_t sampler_lock = PTHREAD_MUTEX_INITIALIZER;
+static pthread_cond_t sampler_cv;
+static bool sampler_stop;
+
+static void *stats_sampler(void *arg)
+{
+	struct timespec next;
+
+	thr_stats_off = true;	/* its own formatting is not the program's */
+	clock_gettime(CLOCK_MONOTONIC, &next);
+	pthread_mutex_lock(&sampler_lock);
+	while (!sampler_stop) {
+		next.tv_nsec += (long)(stats_interval_ms % 1000) * 1000000;
+		next.tv_sec += stats_interval_ms / 1000 + next.tv_nsec / 1000000000;
+		next.tv_nsec %= 1000000000;
+		while (!sampler_stop &&
+		       pthread_cond_timedwait(&sampler_cv, &sampler_lock, &next) != ETIMEDOUT)
+			;
+		if (sampler_stop)
+			break;
+		pthread_mutex_unlock(&sampler_lock);
+		pthread_mutex_lock(&stats_out_lock);
+		stats_emit_interval();
+		pthread_mutex_unlock(&stats_out_lock);
+		pthread_mutex_lock(&sampler_lock);
+	}
+	pthread_mutex_unlock(&sampler_lock);
+	return NULL;
+}
+
+/* Started by the first counted call (not from the constructor or the fork
+ * handler, where creating a thread is not safe), at most once per process. */
+static void stats_start_sampler(void)
+{
+	unsigned char expected = 0;
+	pthread_condattr_t attr;
+	sigset_t all, old;
+
+	if (!atomic_compare_exchange_strong(&sampler_started, &expected, 1))
+		return;
+	pthread_condattr_init(&attr);
+	pthread_condattr_setclock(&attr, CLOCK_MONOTONIC);
+	pthread_cond_init(&sampler_cv, &attr);
+	pthread_condattr_destroy(&attr);
+	sampler_stop = false;
+	/* the sampler must not take the program's signals */
+	sigfillset(&all);
+	pthread_sigmask(SIG_SETMASK, &all, &old);
+	if (pthread_create(&sampler_thread, NULL, stats_sampler, NULL))
+		atomic_store(&sampler_started, 0);
+	pthread_sigmask(SIG_SETMASK, &old, NULL);
+}
+
+static void stats_stop_sampler(void)
+{
+	if (!atomic_load(&sampler_started))
+		return;
+	pthread_mutex_lock(&sampler_lock);
+	sampler_stop = true;
+	pthread_cond_signal(&sampler_cv);
+	pthread_mutex_unlock(&sampler_lock);
+	pthread_join(sampler_thread, NULL);
+	atomic_store(&sampler_started, 0);
+	/* the part of the last interval since the sampler's final record */
+	pthread_mutex_lock(&stats_out_lock);
+	stats_emit_interval();
+	pthread_mutex_unlock(&stats_out_lock);
+}
+
+/* ---- fork ----
+ *
+ * The forking thread holds both stats locks across fork() so the child never
+ * inherits one held by a thread that does not exist there. */
+static void stats_fork_prepare(void)
+{
+	pthread_mutex_lock(&stats_out_lock);
+	pthread_mutex_lock(&stats_lock);
+	pthread_mutex_lock(&sampler_lock);
+}
+
+static void stats_fork_parent(void)
+{
+	pthread_mutex_unlock(&sampler_lock);
+	pthread_mutex_unlock(&stats_lock);
+	pthread_mutex_unlock(&stats_out_lock);
+}
+
+/* In the child: drop every block (their threads are gone, and the child
+ * counts from zero) and forget the parent's sampler. */
+static void stats_fork_child(void)
+{
+	struct dto_thread_stats *t = stats_threads, *next;
+
+	pthread_mutex_unlock(&sampler_lock);
+	pthread_mutex_unlock(&stats_lock);
+	pthread_mutex_unlock(&stats_out_lock);
+	for (; t; t = next) {
+		next = t->next;
+		munmap(t, sizeof(*t));
+	}
+	stats_threads = NULL;
+	orig_memset(&stats_retired, 0, sizeof(stats_retired));
+	orig_memset(&sampler_prev, 0, sizeof(sampler_prev));
+	sampler_prev_ns = 0;
+	if (stats_key_created)
+		pthread_setspecific(stats_key, NULL);
+	thr_stats = NULL;
+	thr_stats_off = false;
+	sampler_stop = false;
+	atomic_store(&sampler_started, 0);
+	stats_out_len = 0;
+	/* with %p the child writes a file of its own: init_dto opens it */
+	if (stats_file_per_pid && stats_fd != -1) {
+		close(stats_fd);
+		stats_fd = -1;
 	}
 }
 #endif
@@ -1454,28 +2035,49 @@ static int init_dto(void)
 		}
 
 		if (collect_stats) {
+			env_str = getenv("DTO_STATS_FORMAT");
+			if (env_str != NULL && !strcmp(env_str, "json"))
+				stats_format = STATS_FORMAT_JSON;
+			else if (env_str != NULL && strcmp(env_str, "text"))
+				LOG_ERROR("Unknown DTO_STATS_FORMAT '%s', using text\n", env_str);
+
+			env_str = getenv("DTO_STATS_INTERVAL_MS");
+			if (env_str != NULL) {
+				char *end;
+
+				errno = 0;
+				stats_interval_ms = strtoul(env_str, &end, 10);
+				/* a day: beyond that it is a typo, and the timespec math stays exact */
+				if (errno || end == env_str || *end || *env_str == '-' ||
+				    stats_interval_ms > 86400000UL) {
+					LOG_ERROR("Invalid DTO_STATS_INTERVAL_MS '%s', no interval records\n", env_str);
+					stats_interval_ms = 0;
+				}
+			}
+
+			/* A forked child keeps the parent's (O_APPEND) file,
+			 * unless the name has a %p (see stats_fork_child). */
 			env_str = getenv("DTO_STATS_FILE");
-                        if (env_str != NULL) {
-                                struct stat st;
-
-                                if (lstat(env_str, &st) == -1 ||
-                                        (st.st_mode & S_IFMT) == S_IFREG) {
-                                        stats_fd = open(env_str, O_WRONLY | O_CREAT | O_TRUNC | O_CLOEXEC, 0600);
-                                        if (stats_fd == -1)
-                                                LOG_ERROR("Failed to open DTO_STATS_FILE '%s': %s\n", env_str, strerror(errno));
-                                } else {
-                                        LOG_ERROR("DTO_STATS_FILE '%s' must reference a regular file\n", env_str);
-                                }
-                        }
-
-                        clock_gettime(CLOCK_BOOTTIME, &dto_start_time);
-                        /* Change the log level to 'trace' so that the
-                         * stats can be logged
-                         */
-                        log_level = LOG_LEVEL_TRACE;
-                        dto_stats_log("Stats collection enabled, some latencies may be higher due to collection overhead.\n");
-                        install_sigint_handler();
-                }
+			if (env_str != NULL && stats_fd == -1)
+				stats_open_file(env_str);
+
+			if (!stats_key_created &&
+			    pthread_key_create(&stats_key, stats_thread_exit) == 0)
+				stats_key_created = true;
+			if (!stats_key_created) {
+				LOG_ERROR("Stats collection disabled: no thread key\n");
+				collect_stats = 0;
+			}
+		}
+
+		if (collect_stats) {
+			clock_gettime(CLOCK_BOOTTIME, &dto_start_time);
+			if (stats_format == STATS_FORMAT_TEXT) {
+				stats_printf("Stats collection enabled, some latencies may be higher due to collection overhead.\n");
+				stats_flush();
+			}
+			install_sigint_handler();
+		}
 #endif
 
 		/* Register fork handler for the child process */
@@ -1485,7 +2087,11 @@ static int init_dto(void)
 			 * take that risk and disable DSA offload for parent
 			 * as well.
 			 */
+#ifdef DTO_STATS_SUPPORT
+			if (!pthread_atfork(stats_fork_prepare, stats_fork_parent, child))
+#else
 			if (!pthread_atfork(NULL, NULL, child))
+#endif
 				fork_handler_registered = 1;
 			else {
 				LOG_ERROR("Setting fork() handler failed. "
@@ -1634,7 +2240,11 @@ static void cleanup_dto(void)
 		close(wqs[i].wq_fd);
 	}
 #ifdef DTO_STATS_SUPPORT
-        print_stats();
+	if (collect_stats) {
+		stats_shutdown = true;
+		stats_stop_sampler();
+		print_stats();
+	}
         restore_sigint_handler();
 #endif
         if (log_fd != -1)
diff --git a/dtoctl.c b/dtoctl.c
index 7c72444..4aa7c50 100644
--- a/dtoctl.c
+++ b/dtoctl.c
@@ -14,14 +14,21 @@
 
 #include <argp.h>
 #include <ctype.h>
+#include <dirent.h>
 #include <errno.h>
 #include <libgen.h>
 #include <limits.h>
+#include <signal.h>
 #include <stdbool.h>
 #include <stdio.h>
 #include <stdlib.h>
 #include <string.h>
+#include <time.h>
 #include <unistd.h>
+#include <sys/resource.h>
+#include <sys/stat.h>
+#include <sys/utsname.h>
+#include <sys/wait.h>
 
 #ifndef LIBDTO_PATH
 #define LIBDTO_PATH "/usr/lib64/libdto.so"
@@ -54,9 +61,17 @@ enum {
 	OPT_NO_AUTO_ADJUST,
 	OPT_STDC_ONLY,
 	OPT_STATS,
+	OPT_STATS_FORMAT,
+	OPT_STATS_INTERVAL,
 	OPT_DRY_RUN,
+	OPT_STATS_RUN,
+	OPT_REPEAT,
+	OPT_PERF,
 };
 
+/* --repeat limit, to catch typos before hours of runs */
+#define REPEAT_MAX		1000
+
 struct opts {
 	int idx;
 	char *exename;
@@ -77,6 +92,8 @@ struct opts {
 	const char *log_file;
 	const char *log_level;
 	const char *stats_file;
+	const char *stats_format;
+	const char *stats_interval;
 
 	/* flag options (false means "not specified") */
 	bool no_memcpy;
@@ -87,6 +104,12 @@ struct opts {
 	bool no_auto_adjust;
 	bool stdc_only;
 	bool stats;
+
+	/* --stats-run */
+	bool stats_run;
+	const char *stats_dir;	/* NULL: a new directory named for the date */
+	unsigned long repeat;	/* 0: not specified */
+	const char *perf;	/* perf binary; NULL: don't profile */
 };
 
 /* (a part of) output in --help and --version options (generated by argp runtime) */
@@ -188,8 +211,19 @@ static struct argp_option dtoctl_options[] = {
 	{.name = "stats-file",
 	 .key = OPT_STATS_FILE,
 	 .arg = "PATH",
-	 .doc = "Write the stats histogram to PATH instead of standard output; "
-		"implies --stats (sets DTO_STATS_FILE; default: standard output)"},
+	 .doc = "Write the stats to PATH instead of standard output; %p in PATH "
+		"becomes the process ID; implies --stats (sets DTO_STATS_FILE; "
+		"default: standard output)"},
+	{.name = "stats-format",
+	 .key = OPT_STATS_FORMAT,
+	 .arg = "FORMAT",
+	 .doc = "Stats output: text tables or json (JSON Lines); implies --stats "
+		"(sets DTO_STATS_FORMAT; default: text)"},
+	{.name = "stats-interval",
+	 .key = OPT_STATS_INTERVAL,
+	 .arg = "MS",
+	 .doc = "Also write the stats for every MS milliseconds; implies --stats "
+		"(sets DTO_STATS_INTERVAL_MS; default: summary at exit only)"},
 	{.name = "log-file",
 	 .key = OPT_LOG_FILE,
 	 .arg = "PATH",
@@ -201,6 +235,27 @@ static struct argp_option dtoctl_options[] = {
 	 .doc = "Verbosity of DTO logging: 0, 1, or 2 (sets DTO_LOG_LEVEL; "
 		"default: 0)"},
 
+	{.doc = "Stats run (how much DSA speeds the program up):", .group = 6},
+	{.name = "stats-run",
+	 .key = OPT_STATS_RUN,
+	 .arg = "DIR",
+	 .flags = OPTION_ARG_OPTIONAL,
+	 .doc = "Run the program CPU-only, then with DSA, collecting stats as JSON "
+		"in DIR (default: a new dto-stats-<date>-<time>), then summarize the "
+		"speedup with dto-stats-report"},
+	{.name = "repeat",
+	 .key = OPT_REPEAT,
+	 .arg = "N",
+	 .doc = "With --stats-run, run each way N times; the report takes the "
+		"median (default: 1)"},
+	{.name = "perf",
+	 .key = OPT_PERF,
+	 .arg = "PERF",
+	 .flags = OPTION_ARG_OPTIONAL,
+	 .doc = "With --stats-run, also profile one more run each way with "
+		"\"PERF record\" (default: perf from PATH), to show where the time "
+		"goes: copies, memset, memcmp, CRC, kernel copies, DTO itself"},
+
 	{0},
 };
 
@@ -374,6 +429,35 @@ static error_t parse_option(int key, char *arg, struct argp_state *state)
 	case OPT_STATS_FILE:
 		opts->stats_file = check_path(state, "stats-file", arg);
 		break;
+	case OPT_STATS_FORMAT:
+		if (strcmp(arg, "text") && strcmp(arg, "json"))
+			argp_error(state, "invalid stats format '%s': must be text or json",
+				   arg);
+		opts->stats_format = arg;
+		break;
+	case OPT_STATS_INTERVAL:
+		if (!parse_ulong(arg, &val) || val == 0 || val > INT_MAX)
+			argp_error(state, "invalid stats interval '%s': must be a positive "
+				   "number of milliseconds", arg);
+		opts->stats_interval = arg;
+		break;
+	case OPT_STATS_RUN:
+		opts->stats_run = true;
+		opts->stats_dir = arg ? check_path(state, "stats-run", arg) : NULL;
+		/* libdto would expand a %p in the stats file names */
+		if (arg && strstr(arg, "%p"))
+			argp_error(state, "--stats-run directory '%s' cannot contain '%%p'",
+				   arg);
+		break;
+	case OPT_REPEAT:
+		if (!parse_ulong(arg, &val) || val == 0 || val > REPEAT_MAX)
+			argp_error(state, "invalid repeat count '%s': must be from 1 to %d",
+				   arg, REPEAT_MAX);
+		opts->repeat = val;
+		break;
+	case OPT_PERF:
+		opts->perf = arg ? check_path(state, "perf", arg) : "perf";
+		break;
 
 	case ARGP_KEY_ARG:
 		/* leave the program and its arguments to ARGP_KEY_ARGS */
@@ -387,8 +471,17 @@ static error_t parse_option(int key, char *arg, struct argp_state *state)
 
 	case ARGP_KEY_END:
 		/* --dry-run can show the environment without a program */
-		if (opts->exename == NULL && !opts->dry_run)
+		if (opts->exename == NULL && (!opts->dry_run || opts->stats_run))
 			argp_usage(state);
+		if ((opts->repeat || opts->perf) && !opts->stats_run)
+			argp_error(state, "--repeat and --perf need --stats-run");
+		/* a stats run chooses these itself, for each run */
+		if (opts->stats_run && opts->stdc_only)
+			argp_error(state, "--stats-run runs the program CPU-only itself; "
+				   "drop --stdc-only");
+		if (opts->stats_run && (opts->stats_file || opts->stats_format))
+			argp_error(state, "--stats-run writes JSON stats to its own "
+				   "directory; drop --stats-file and --stats-format");
 		break;
 
 	default:
@@ -467,6 +560,12 @@ static void setup_dto_environ(struct opts *opts)
 		set_option_var("DTO_LOG_LEVEL", opts->log_level);
 	if (opts->stats_file)
 		set_option_var("DTO_STATS_FILE", opts->stats_file);
+	if (opts->stats_format)
+		set_option_var("DTO_STATS_FORMAT", opts->stats_format);
+	if (opts->stats_run)
+		set_option_var("DTO_STATS_FORMAT", "json");
+	if (opts->stats_interval)
+		set_option_var("DTO_STATS_INTERVAL_MS", opts->stats_interval);
 
 	/* flag options: only set when explicitly specified */
 	if (opts->no_memcpy)
@@ -483,8 +582,9 @@ static void setup_dto_environ(struct opts *opts)
 		set_option_var("DTO_AUTO_ADJUST_KNOBS", "0");
 	if (opts->stdc_only)
 		set_option_var("DTO_USESTDC_CALLS", "1");
-	/* libdto ignores DTO_STATS_FILE unless stats are collected */
-	if (opts->stats || opts->stats_file)
+	/* libdto ignores the other DTO_STATS_* unless stats are collected */
+	if (opts->stats || opts->stats_file || opts->stats_format ||
+	    opts->stats_interval || opts->stats_run)
 		set_option_var("DTO_COLLECT_STATS", "1");
 }
 
@@ -673,6 +773,479 @@ static void print_dry_run(const char *lib, const char *source, char **argv)
 	putchar('\n');
 }
 
+/*
+ * --stats-run: the program runs with DTO_USESTDC_CALLS set to each of these
+ * in turn, with JSON stats going to <dir>/<run>.jsonl.  Runs are named
+ * <mode>-<round>, or <mode>-perf for the profiled ones.
+ */
+struct run_mode {
+	const char *name;
+	const char *usestdc;
+};
+
+static const struct run_mode run_modes[] = {
+	{ "cpu", "1" },
+	{ "dsa", "0" },
+};
+
+#define NR_RUN_MODES	(sizeof(run_modes) / sizeof(run_modes[0]))
+
+/*
+ * perf record -q -e cpu-clock -F 999 -o <data> -- <dtoctl> -l <lib> -- <program>
+ * cpu-clock samples time on CPU, which is what the report divides up, and
+ * needs no hardware counters (so it also works in VMs).
+ */
+static const char *const perf_record[] = {
+	"record", "-q", "-e", "cpu-clock", "-F", "999", "-o",
+};
+
+#define PERF_RECORD_ARGS	(sizeof(perf_record) / sizeof(perf_record[0]))
+#define PERF_ARGV_OUTPUT	(1 + PERF_RECORD_ARGS)
+
+/* LD_PRELOAD as dtoctl found it, for perf and dto-stats-report */
+static char *orig_preload;
+
+/* The caller's LD_PRELOAD without libdto, for perf and the report: they
+ * must not be measured, and setup_preload may have kept a libdto the
+ * caller had exported. */
+static void restore_preload(void)
+{
+	char *copy = NULL, *preload = NULL, *entry, *save;
+	int err;
+
+	if (orig_preload == NULL) {
+		err = unsetenv("LD_PRELOAD");
+	} else {
+		copy = strdup(orig_preload);
+		preload = malloc(strlen(orig_preload) + 1);
+		if (copy == NULL || preload == NULL) {
+			perror("dtoctl: cannot build LD_PRELOAD");
+			exit(EXIT_DTOCTL_FAILURE);
+		}
+		*preload = '\0';
+		for (entry = strtok_r(copy, ": ", &save); entry;
+		     entry = strtok_r(NULL, ": ", &save)) {
+			if (is_libdto(entry))
+				continue;
+			if (*preload)
+				strcat(preload, ":");
+			strcat(preload, entry);
+		}
+		err = *preload ? setenv("LD_PRELOAD", preload, 1) :
+				 unsetenv("LD_PRELOAD");
+	}
+	free(copy);
+	free(preload);
+	if (err) {
+		perror("dtoctl: LD_PRELOAD");
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+}
+
+static char *self_exe(void)
+{
+	static char exe[PATH_MAX];
+	ssize_t len;
+
+	len = readlink("/proc/self/exe", exe, sizeof(exe) - 1);
+	if (len < 0) {
+		perror("dtoctl: /proc/self/exe");
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	exe[len] = '\0';
+	return exe;
+}
+
+/* true if cmd names a program execvp() would run */
+static bool executable(const char *cmd)
+{
+	const char *path = getenv("PATH");
+	char *dirs, *dir, *save, *file;
+	bool found = false;
+
+	if (strchr(cmd, '/'))
+		return access(cmd, X_OK) == 0;
+	dirs = strdup(path ? path : "/bin:/usr/bin");
+	for (dir = strtok_r(dirs, ":", &save); dir && !found;
+	     dir = strtok_r(NULL, ":", &save)) {
+		if (asprintf(&file, "%s/%s", *dir ? dir : ".", cmd) < 0)
+			break;
+		found = access(file, X_OK) == 0;
+		free(file);
+	}
+	free(dirs);
+	return found;
+}
+
+/* The command for a profiled run; the output file is filled in per run. */
+static char **perf_command(const char *perf, const char *lib, char **argv)
+{
+	int argc = 0, i = 0;
+	char **cmd;
+
+	while (argv[argc])
+		argc++;
+	cmd = calloc(PERF_RECORD_ARGS + argc + 8, sizeof(*cmd));
+	if (cmd == NULL) {
+		perror("dtoctl");
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	cmd[i++] = (char *)perf;
+	for (size_t j = 0; j < PERF_RECORD_ARGS; j++)
+		cmd[i++] = (char *)perf_record[j];
+	cmd[i++] = NULL;	/* PERF_ARGV_OUTPUT */
+	cmd[i++] = "--";
+	/* run dtoctl again under perf, so DTO is preloaded into the program only */
+	cmd[i++] = self_exe();
+	cmd[i++] = "-l";
+	cmd[i++] = (char *)lib;
+	cmd[i++] = "--";
+	memcpy(cmd + i, argv, argc * sizeof(*cmd));
+	return cmd;
+}
+
+static void run_name(char *buf, size_t size, const struct run_mode *mode,
+		     unsigned long round)
+{
+	if (round)
+		snprintf(buf, size, "%s-%lu", mode->name, round);
+	else
+		snprintf(buf, size, "%s-perf", mode->name);
+}
+
+static void json_string(FILE *f, const char *s)
+{
+	fputc('"', f);
+	for (; *s; s++) {
+		unsigned char c = *s;
+
+		if (c == '"' || c == '\\')
+			fprintf(f, "\\%c", c);
+		else if (c < 0x20)
+			fprintf(f, "\\u%04x", c);
+		else
+			fputc(c, f);
+	}
+	fputc('"', f);
+}
+
+/* the variables each run sets for itself */
+static bool is_run_var(const char *name, int len)
+{
+	return (len == 17 && !strncmp(name, "DTO_USESTDC_CALLS", len)) ||
+	       (len == 14 && !strncmp(name, "DTO_STATS_FILE", len));
+}
+
+/* First line of runs.jsonl: what was run and how. */
+static void write_run_header(FILE *runs, const struct opts *opts, const char *lib,
+			     unsigned long repeat, char **argv)
+{
+	time_t now = time(NULL);
+	struct utsname un;
+	char start[32], *name;
+	bool first = true;
+	char **env;
+	int len;
+
+	strftime(start, sizeof(start), "%Y-%m-%dT%H:%M:%S%z", localtime(&now));
+	if (uname(&un))
+		strcpy(un.nodename, "");
+
+	fputs("{\"type\":\"stats_run\",\"dtoctl_version\":", runs);
+	json_string(runs, DTO_VERSION);
+	fputs(",\"start\":", runs);
+	json_string(runs, start);
+	fputs(",\"host\":", runs);
+	json_string(runs, un.nodename);
+	fputs(",\"library\":", runs);
+	json_string(runs, lib);
+	fprintf(runs, ",\"repeat\":%lu,\"perf\":", repeat);
+	if (opts->perf)
+		json_string(runs, opts->perf);
+	else
+		fputs("null", runs);
+	fputs(",\"command\":[", runs);
+	for (int i = 0; argv[i]; i++) {
+		if (i)
+			fputc(',', runs);
+		json_string(runs, argv[i]);
+	}
+	fputs("],\"env\":{", runs);
+	for (env = environ; *env; env++) {
+		len = strcspn(*env, "=");
+		if (strncmp(*env, "DTO_", 4) || (*env)[len] != '=' ||
+		    is_run_var(*env, len))
+			continue;
+		name = strndup(*env, len);
+		if (name == NULL)
+			continue;
+		if (!first)
+			fputc(',', runs);
+		json_string(runs, name);
+		fputc(':', runs);
+		json_string(runs, *env + len + 1);
+		free(name);
+		first = false;
+	}
+	fputs("}}\n", runs);
+	fflush(runs);
+}
+
+/* Make the --stats-run directory, or check that it is empty; returns its absolute path. */
+static char *make_stats_dir(const char *dir)
+{
+	struct dirent *ent;
+	char *path;
+	DIR *d;
+
+	/* libdto expands %p in the stats file names; the working directory
+	 * may bring one in, so check the full path before creating anything */
+	if (dir[0] != '/' && (path = getcwd(NULL, 0)) != NULL) {
+		if (strstr(path, "%p")) {
+			fprintf(stderr, "dtoctl: stats directory '%s/%s' cannot contain '%%p'\n",
+				path, dir);
+			exit(EXIT_DTOCTL_FAILURE);
+		}
+		free(path);
+	}
+	if (mkdir(dir, 0777) && errno != EEXIST) {
+		fprintf(stderr, "dtoctl: cannot create '%s': %s\n", dir, strerror(errno));
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	/* one run's results only, so the report can't mix two */
+	d = opendir(dir);
+	if (d == NULL) {
+		fprintf(stderr, "dtoctl: cannot open '%s': %s\n", dir, strerror(errno));
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	while ((ent = readdir(d))) {
+		if (strcmp(ent->d_name, ".") && strcmp(ent->d_name, "..")) {
+			fprintf(stderr, "dtoctl: '%s' is not empty; give --stats-run a "
+				"new directory\n", dir);
+			exit(EXIT_DTOCTL_FAILURE);
+		}
+	}
+	closedir(d);
+
+	path = realpath(dir, NULL);
+	if (path == NULL) {
+		fprintf(stderr, "dtoctl: '%s': %s\n", dir, strerror(errno));
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	if (strstr(path, "%p")) {
+		fprintf(stderr, "dtoctl: stats directory '%s' cannot contain '%%p'\n", path);
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	return path;
+}
+
+static double seconds(struct timeval tv)
+{
+	return tv.tv_sec + tv.tv_usec / 1e6;
+}
+
+/*
+ * One run of the program for --stats-run, under perf if perf_argv is set.
+ * Appends what happened to runs.jsonl.  Returns 0 if the program succeeded,
+ * else the status for dtoctl to exit with.
+ */
+static int stats_run_one(FILE *runs, const char *dir, const struct run_mode *mode,
+			 unsigned long round, char **perf_argv, char **argv)
+{
+	struct sigaction ign = { .sa_handler = SIG_IGN }, old_int, old_quit;
+	char name[32], *stats, *data = NULL;
+	struct timespec start, end;
+	sigset_t block, saved;
+	struct rusage ru;
+	int status, err, ret;
+	pid_t pid;
+
+	run_name(name, sizeof(name), mode, round);
+	if (asprintf(&stats, "%s/%s.jsonl", dir, name) < 0 ||
+	    (perf_argv && asprintf(&data, "%s/%s.data", dir, name) < 0)) {
+		perror("dtoctl");
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	if (perf_argv) {
+		perf_argv[PERF_ARGV_OUTPUT] = data;
+		argv = perf_argv;
+	}
+	fprintf(stderr, "dtoctl: run %s\n", name);
+
+	/* as system() does: Ctrl-C and Ctrl-\ stop the program, not dtoctl */
+	sigemptyset(&block);
+	sigaddset(&block, SIGINT);
+	sigaddset(&block, SIGQUIT);
+	sigprocmask(SIG_BLOCK, &block, &saved);
+
+	clock_gettime(CLOCK_MONOTONIC, &start);
+	pid = fork();
+	if (pid == 0) {
+		sigprocmask(SIG_SETMASK, &saved, NULL);
+		if (setenv("DTO_USESTDC_CALLS", mode->usestdc, 1) ||
+		    setenv("DTO_STATS_FILE", stats, 1))
+			_exit(EXIT_DTOCTL_FAILURE);
+		/* under perf, the inner dtoctl sets LD_PRELOAD for the program */
+		if (perf_argv)
+			restore_preload();
+		execvp(argv[0], argv);
+		err = errno;
+		fprintf(stderr, "dtoctl: cannot run '%s': %s\n", argv[0], strerror(err));
+		_exit(err == ENOENT ? EXIT_ENOENT : EXIT_CANNOT_INVOKE);
+	}
+	if (pid < 0) {
+		perror("dtoctl: fork");
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	sigaction(SIGINT, &ign, &old_int);
+	sigaction(SIGQUIT, &ign, &old_quit);
+	sigprocmask(SIG_SETMASK, &saved, NULL);
+
+	while (wait4(pid, &status, 0, &ru) < 0) {
+		if (errno != EINTR) {
+			perror("dtoctl: wait4");
+			exit(EXIT_DTOCTL_FAILURE);
+		}
+	}
+	clock_gettime(CLOCK_MONOTONIC, &end);
+	sigaction(SIGINT, &old_int, NULL);
+	sigaction(SIGQUIT, &old_quit, NULL);
+
+	fprintf(runs, "{\"type\":\"run\",\"name\":\"%s\",\"mode\":\"%s\",\"round\":%lu,"
+		"\"profiled\":%s,\"stats\":\"%s.jsonl\",", name, mode->name, round,
+		perf_argv ? "true" : "false", name);
+	if (perf_argv)
+		fprintf(runs, "\"perf_data\":\"%s.data\",", name);
+	fprintf(runs, "\"wall_s\":%.6f,\"user_s\":%.6f,\"sys_s\":%.6f,\"maxrss_kb\":%ld,",
+		(end.tv_sec - start.tv_sec) + (end.tv_nsec - start.tv_nsec) / 1e9,
+		seconds(ru.ru_utime), seconds(ru.ru_stime), ru.ru_maxrss);
+	if (WIFEXITED(status)) {
+		ret = WEXITSTATUS(status);
+		fprintf(runs, "\"exit_status\":%d}\n", ret);
+		if (ret)
+			fprintf(stderr, "dtoctl: run %s exited with status %d%s\n",
+				name, ret, perf_argv ? "" : "; stopping");
+	} else {
+		ret = 128 + WTERMSIG(status);
+		fprintf(runs, "\"signal\":%d}\n", WTERMSIG(status));
+		fprintf(stderr, "dtoctl: run %s was killed by signal %d%s\n",
+			name, WTERMSIG(status), perf_argv ? "" : "; stopping");
+	}
+	fflush(runs);
+	free(stats);
+	free(data);
+	return ret;
+}
+
+/* Summarize the runs with dto-stats-report, preferring the one beside dtoctl. */
+static int run_report(const char *dir)
+{
+	char *report;
+
+	restore_preload();
+	if (asprintf(&report, "%s/dto-stats-report", dirname(self_exe())) >= 0 &&
+	    access(report, X_OK) == 0)
+		execl(report, report, dir, (char *)NULL);
+	execlp("dto-stats-report", "dto-stats-report", dir, (char *)NULL);
+	fprintf(stderr, "dtoctl: stats saved in %s; dto-stats-report (not found) "
+		"summarizes them\n", dir);
+	return 0;
+}
+
+static void print_run(const char *dir, const struct run_mode *mode,
+		      unsigned long round, char **argv)
+{
+	char name[32], *stats;
+
+	run_name(name, sizeof(name), mode, round);
+	if (asprintf(&stats, "%s/%s.jsonl", dir, name) < 0)
+		return;
+	printf("run %s: DTO_USESTDC_CALLS=%s DTO_STATS_FILE=", name, mode->usestdc);
+	print_quoted(stats);
+	for (int i = 0; argv[i]; i++) {
+		putchar(' ');
+		print_quoted(argv[i]);
+	}
+	putchar('\n');
+	free(stats);
+}
+
+/* --stats-run: run the program CPU-only and with DSA, then report. */
+static int stats_run(struct opts *opts, const char *lib, const char *source,
+		     char **argv)
+{
+	unsigned long repeat = opts->repeat ? opts->repeat : 1;
+	char **perf_argv = NULL, *dir, *data;
+	time_t now = time(NULL);
+	char name[64];
+	FILE *runs;
+	size_t m;
+	int ret;
+
+	if (opts->stats_dir) {
+		dir = (char *)opts->stats_dir;
+	} else {
+		strftime(name, sizeof(name), "dto-stats-%Y%m%d-%H%M%S", localtime(&now));
+		dir = name;
+	}
+	if (opts->perf)
+		perf_argv = perf_command(opts->perf, lib, argv);
+
+	if (opts->dry_run) {
+		print_dry_run(lib, source, NULL);
+		printf("stats run: %s\n", dir);
+		for (unsigned long round = 1; round <= repeat; round++)
+			for (m = 0; m < NR_RUN_MODES; m++)
+				print_run(dir, &run_modes[m], round, argv);
+		for (m = 0; perf_argv && m < NR_RUN_MODES; m++) {
+			if (asprintf(&data, "%s/%s-perf.data", dir, run_modes[m].name) < 0)
+				break;
+			perf_argv[PERF_ARGV_OUTPUT] = data;
+			print_run(dir, &run_modes[m], 0, perf_argv);
+			free(data);
+		}
+		printf("report: dto-stats-report ");
+		print_quoted(dir);
+		putchar('\n');
+		return 0;
+	}
+
+	/* fail now rather than after the timed runs */
+	if (opts->perf && !executable(opts->perf)) {
+		fprintf(stderr, "dtoctl: cannot find perf '%s'\n", opts->perf);
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	dir = make_stats_dir(dir);
+	if (asprintf(&data, "%s/runs.jsonl", dir) < 0 ||
+	    (runs = fopen(data, "we")) == NULL) {
+		fprintf(stderr, "dtoctl: cannot write runs.jsonl in '%s': %s\n", dir,
+			strerror(errno));
+		exit(EXIT_DTOCTL_FAILURE);
+	}
+	free(data);
+	write_run_header(runs, opts, lib, repeat, argv);
+
+	/* alternate the modes, so drift over time affects both alike */
+	for (unsigned long round = 1; round <= repeat; round++) {
+		for (m = 0; m < NR_RUN_MODES; m++) {
+			ret = stats_run_one(runs, dir, &run_modes[m], round, NULL, argv);
+			if (ret)
+				return ret;
+		}
+	}
+	/* the profiles are extras: the timed runs are reported without them */
+	for (m = 0; perf_argv && m < NR_RUN_MODES; m++) {
+		ret = stats_run_one(runs, dir, &run_modes[m], 0, perf_argv, argv);
+		if (ret) {
+			fprintf(stderr, "dtoctl: reporting the timed runs without a profile\n");
+			break;
+		}
+	}
+	fclose(runs);
+	fflush(stdout);
+	return run_report(dir);
+}
+
 int main(int argc, char *argv[])
 {
 	struct opts opts = {0};
@@ -684,7 +1257,9 @@ int main(int argc, char *argv[])
 		       "\vOnly the options given are applied: inherited DTO_* variables "
 		       "are kept unless an option overrides them or -i is given.  "
 		       "Exit status is the program's own, or 125 if dtoctl fails, "
-		       "126 if the program cannot be run, 127 if it is not found.",
+		       "126 if the program cannot be run, 127 if it is not found.  "
+		       "With --stats-run, a run that fails stops the rest, and "
+		       "dtoctl exits with its status.",
 	};
 	const char *source;
 	char *lib;
@@ -698,8 +1273,13 @@ int main(int argc, char *argv[])
 	setup_dto_environ(&opts);
 
 	lib = resolve_library(&opts, &source);
+	if (opts.stats_run && getenv("LD_PRELOAD"))
+		orig_preload = strdup(getenv("LD_PRELOAD"));
 	setup_preload(lib);
 
+	if (opts.stats_run)
+		return stats_run(&opts, lib, source, argv + opts.idx);
+
 	if (opts.dry_run) {
 		print_dry_run(lib, source, opts.exename ? argv + opts.idx : NULL);
 		return 0;
diff --git a/tests/test_dtoctl.sh b/tests/test_dtoctl.sh
index 8cf4cf0..f9d3a26 100755
--- a/tests/test_dtoctl.sh
+++ b/tests/test_dtoctl.sh
@@ -54,14 +54,16 @@ status() {
 run "$DTOCTL" -l "$LIB" --dry-run -w umwait -b 64K -c 0.25 -n 2 -q 'wq0.0;wq2.0' \
     --umwait-delay 5000 --overlapping-memmove=dsa --no-memcpy --no-memmove \
     --no-memset --no-memcmp --no-cache-control --no-auto-adjust --stdc-only \
-    --log-file=/tmp/dto.log --log-level 1 --stats-file=/tmp/dto.stats prog 'a b'
+    --log-file=/tmp/dto.log --log-level 1 --stats-file=/tmp/dto.stats \
+    --stats-format json --stats-interval 10 prog 'a b'
 status "dry-run succeeds" 0
 for v in DTO_WAIT_METHOD=umwait DTO_MIN_BYTES=65536 DTO_CPU_SIZE_FRACTION=0.25 \
          DTO_IS_NUMA_AWARE=2 "DTO_WQ_LIST='wq0.0;wq2.0'" DTO_UMWAIT_DELAY=5000 \
          DTO_OVERLAPPING_MEMMOVE_ACTION=1 DTO_DSA_MEMCPY=0 DTO_DSA_MEMMOVE=0 \
          DTO_DSA_MEMSET=0 DTO_DSA_MEMCMP=0 DTO_DSA_CC=0 DTO_AUTO_ADJUST_KNOBS=0 \
          DTO_USESTDC_CALLS=1 DTO_LOG_FILE=/tmp/dto.log DTO_LOG_LEVEL=1 \
-         DTO_STATS_FILE=/tmp/dto.stats DTO_COLLECT_STATS=1; do
+         DTO_STATS_FILE=/tmp/dto.stats DTO_STATS_FORMAT=json DTO_STATS_INTERVAL_MS=10 \
+         DTO_COLLECT_STATS=1; do
     has "option sets $v" "^$v  # option\$"
 done
 has "dry-run shows the library" "^library: $LIB (--library)\$"
@@ -111,7 +113,10 @@ has "another libdto in LD_PRELOAD is dropped" "dropping '$TMP/old/libdto.so.1'"
 
 # --- invalid values are rejected ---
 for bad in "-w bogus" "-b 64Q" "-b abc" "-b -1" "-c 1.5" "-c nan" "-n 3" \
-           "--log-level x" "--umwait-delay 0" "--overlapping-memmove=gpu"; do
+           "--log-level x" "--umwait-delay 0" "--overlapping-memmove=gpu" \
+           "--stats-format xml" "--stats-interval 0" "--repeat 2" "--perf" \
+           "--stats-run --stdc-only" "--stats-run --stats-file=x" "--stats-run=a%pb" \
+           "--stats-run --repeat 0"; do
     run "$DTOCTL" -l "$LIB" $bad true
     status "rejects $bad" 64
 done
@@ -139,6 +144,56 @@ status "program's exit status is returned" 7
 run "$DTOCTL" --version
 has "--version" "^dtoctl (DTO) "
 
+# --- stats run ---
+run "$DTOCTL" -l "$LIB" --dry-run --stats-run=out --repeat 2 --perf=/usr/bin/perf prog 'a b'
+status "stats run dry-run succeeds" 0
+has "stats run collects stats" "^DTO_COLLECT_STATS=1  # option\$"
+has "stats run collects JSON" "^DTO_STATS_FORMAT=json  # option\$"
+has "first run is CPU-only" \
+    "^run cpu-1: DTO_USESTDC_CALLS=1 DTO_STATS_FILE=out/cpu-1.jsonl prog 'a b'\$"
+has "second run uses DSA" \
+    "^run dsa-1: DTO_USESTDC_CALLS=0 DTO_STATS_FILE=out/dsa-1.jsonl prog 'a b'\$"
+has "--repeat 2 runs each way twice" "^run dsa-2: "
+lacks "--repeat 2 stops at two" "^run cpu-3: "
+has "profiled run goes through perf and dtoctl" \
+    "^run cpu-perf: DTO_USESTDC_CALLS=1 DTO_STATS_FILE=out/cpu-perf.jsonl /usr/bin/perf record .* -o out/cpu-perf.data -- [^ ]*dtoctl -l $LIB -- prog 'a b'\$"
+has "the report comes last" "^report: dto-stats-report out\$"
+
+run "$DTOCTL" -l "$LIB" --dry-run --stats-run
+status "--stats-run needs a program, even with --dry-run" 64
+
+# dash leaves with _exit(), skipping DTO's summary, so end in a program that exits
+run "$DTOCTL" -l "$LIB" --stats-run="$TMP/run" --repeat 2 \
+    sh -c 'echo "env $DTO_USESTDC_CALLS ${DTO_STATS_FILE##*/} $DTO_STATS_FORMAT"; exec true'
+if command -v python3 >/dev/null; then
+    status "stats run succeeds" 0
+else
+    status "stats run succeeds (the report needs python3)" 127
+fi
+has "CPU-only runs set DTO_USESTDC_CALLS=1" "^env 1 cpu-2.jsonl json\$"
+has "DSA runs set DTO_USESTDC_CALLS=0" "^env 0 dsa-2.jsonl json\$"
+if [ "$(grep -c '"type":"run".*"exit_status":0' "$TMP/run/runs.jsonl")" -eq 4 ] &&
+   grep -q '"type":"summary"' "$TMP/run/dsa-2.jsonl"; then
+    echo "PASS: runs.jsonl records each run, and each run's stats are saved"
+else
+    fail "runs.jsonl records each run, and each run's stats are saved"
+fi
+if command -v python3 >/dev/null; then
+    has "dto-stats-report summarizes the runs" "^DTO stats run: $TMP/run\$"
+fi
+
+run "$DTOCTL" -l "$LIB" --stats-run="$TMP/run" true
+status "a directory with results in it is refused" 125
+
+run "$DTOCTL" -l "$LIB" --stats-run="$TMP/fail" sh -c 'exit 3'
+status "a failing run ends the stats run with its status" 3
+has "the failing run is named" "run cpu-1 exited with status 3; stopping"
+if [ -e "$TMP/fail/dsa-1.jsonl" ]; then
+    fail "no runs after a failure"
+else
+    echo "PASS: no runs after a failure"
+fi
+
 # --- libdto applies the settings (DSA is optional) ---
 run "$DTOCTL" -l "$LIB" --log-level 2 -w yield -b 8K -c 0.5 --no-auto-adjust true
 has "libdto applies the options" \
diff --git a/tests/test_stats.py b/tests/test_stats.py
new file mode 100644
index 0000000..0043159
--- /dev/null
+++ b/tests/test_stats.py
@@ -0,0 +1,148 @@
+#!/usr/bin/env python3
+# SPDX-License-Identifier: MIT
+#
+# Checks DTO_COLLECT_STATS output against test_stats_workload's known calls.
+# Needs no DSA hardware: without it every call is counted on the CPU path.
+#
+# Usage: test_stats.py <test_stats_workload> <libdto.so>
+
+import json
+import os
+import subprocess
+import sys
+import tempfile
+
+WORKLOAD = sys.argv[1]
+LIBDTO = os.path.abspath(sys.argv[2])
+SIZES = {"memcpy": 131072 + 123, "memset": 65536 + 4096,
+         "memmove": 98304 + 7, "memcmp": 262144 + 64}
+SIZE_CHILD = 524288 + 1
+failures = 0
+
+
+def check(desc, ok, detail=""):
+    global failures
+    print(("PASS: " if ok else "FAIL: ") + desc + ("" if ok else "  " + detail))
+    failures += not ok
+
+
+def run(args, env_extra, stats_path):
+    env = dict(os.environ, DTO_COLLECT_STATS="1", DTO_STATS_FILE=stats_path,
+               **env_extra)
+    out = subprocess.run([WORKLOAD] + args, env=env, check=True,
+                         capture_output=True, text=True).stdout
+    with open(stats_path) as f:
+        return out, f.read()
+
+
+def records(text):
+    return [json.loads(line) for line in text.splitlines() if line.strip()]
+
+
+def op_sum(ops, op, field):
+    return sum(g[field] for g in ops.get(op, {}).values())
+
+
+def bucket_for(summary, size):
+    for b in summary["histogram"]:
+        if b["min"] <= size and ("max" not in b or size <= b["max"]):
+            return b
+    return None
+
+
+tmp = tempfile.mkdtemp()
+
+# --- totals and histogram, JSON ---
+threads, iters = 4, 50
+_, text = run(["basic", str(threads), str(iters)], {"DTO_STATS_FORMAT": "json"},
+              os.path.join(tmp, "basic.json"))
+recs = records(text)
+summaries = [r for r in recs if r["type"] == "summary"]
+check("one summary record", len(summaries) == 1, str(len(summaries)))
+s = summaries[0]
+n = threads * iters
+for op, size in SIZES.items():
+    b = bucket_for(s, size)
+    check(f"{op} lands in the bucket for {size} bytes", b is not None and op in b["ops"],
+          json.dumps(b))
+    if b is None or op not in b["ops"]:
+        continue
+    calls = op_sum(b["ops"], op, "calls")
+    failed = b["ops"][op].get("dsa_failed", {}).get("calls", 0)
+    # a call the device left partly done is counted as dsa_failed plus cpu
+    check(f"{op} calls ({n})", calls - failed == n, f"calls {calls} failed {failed}")
+    check(f"{op} bytes ({n * size})", op_sum(b["ops"], op, "bytes") == n * size,
+          str(op_sum(b["ops"], op, "bytes")))
+    check(f"{op} time recorded", op_sum(b["ops"], op, "ns") > 0)
+check("summary carries the configuration", "wait_method" in s["config"])
+
+# --- interval records add up to the totals ---
+_, text = run(["basic", "2", "40", "2"],
+              {"DTO_STATS_FORMAT": "json", "DTO_STATS_INTERVAL_MS": "10"},
+              os.path.join(tmp, "interval.json"))
+recs = records(text)
+intervals = [r for r in recs if r["type"] == "interval"]
+s = [r for r in recs if r["type"] == "summary"][0]
+check("interval records emitted", len(intervals) >= 5, str(len(intervals)))
+check("interval times increase",
+      all(b["t_ms"] > a["t_ms"] for a, b in zip(intervals, intervals[1:])))
+for op in SIZES:
+    for field in ("calls", "bytes"):
+        total = sum(op_sum(r["ops"], op, field) for r in intervals)
+        check(f"intervals sum to the {op} {field} total",
+              total == op_sum(s["totals"], op, field),
+              f"{total} vs {op_sum(s['totals'], op, field)}")
+
+# --- short-lived threads: all counted, blocks released ---
+out, text = run(["churn", "200", "5"], {"DTO_STATS_FORMAT": "json"},
+                os.path.join(tmp, "churn.json"))
+s = [r for r in records(text) if r["type"] == "summary"][0]
+b = bucket_for(s, SIZES["memcpy"])
+calls = op_sum(b["ops"], "memcpy", "calls") - \
+    b["ops"]["memcpy"].get("dsa_failed", {}).get("calls", 0)
+check("calls from 200 exited threads counted", calls == 200 * 5, str(calls))
+before, after = (int(x) for x in out.split()[1::2])
+check("exited threads' stats blocks are unmapped", after - before < 20,
+      f"mappings {before} -> {after}")
+
+# --- fork: separate records per process, appended to one file ---
+out, text = run(["fork", "30"], {"DTO_STATS_FORMAT": "json"},
+                os.path.join(tmp, "fork.json"))
+child = int(out.split()[1])
+by_pid = {r["pid"]: r for r in records(text) if r["type"] == "summary"}
+check("parent and child both report, in one file", len(by_pid) == 2, str(list(by_pid)))
+if child in by_pid and len(by_pid) == 2:
+    parent = [p for p in by_pid if p != child][0]
+    cb = bucket_for(by_pid[child], SIZE_CHILD)
+    check("child counts only its own copies",
+          cb is not None and bucket_for(by_pid[child], SIZES["memcpy"]) is None)
+    pb = bucket_for(by_pid[parent], SIZES["memcpy"])
+    check("parent counts its copies from before and after the fork",
+          pb is not None and op_sum(pb["ops"], "memcpy", "bytes") == 60 * SIZES["memcpy"])
+
+# --- %p: a file per process ---
+pattern = os.path.join(tmp, "pp.%p.json")
+env = dict(os.environ, DTO_COLLECT_STATS="1", DTO_STATS_FILE=pattern,
+           DTO_STATS_FORMAT="json")
+subprocess.run([WORKLOAD, "fork", "5"], env=env, check=True, capture_output=True)
+files = [f for f in os.listdir(tmp) if f.startswith("pp.")]
+check("%p gives each process its own file", len(files) == 2, str(files))
+
+# --- a program that doesn't link libpthread (before glibc 2.34, it isn't libc) ---
+env = dict(os.environ, LD_PRELOAD=LIBDTO, DTO_COLLECT_STATS="1",
+           DTO_STATS_INTERVAL_MS="1", DTO_STATS_FILE=os.path.join(tmp, "true.txt"))
+p = subprocess.run(["true"], env=env, capture_output=True, text=True)
+check("stats work in a program without libpthread", p.returncode == 0, p.stderr.strip())
+
+# --- text format ---
+_, text = run(["basic", "1", "10"], {"DTO_STATS_FORMAT": "text"},
+              os.path.join(tmp, "basic.txt"))
+check("text summary has the operation table", "Number of Memory Operations" in text
+      and "Byte Range" in text)
+check("text summary lists the memcpy bucket", "  131072-163839" in text)
+
+print()
+if failures:
+    print(f"{failures} stats test(s) failed")
+    sys.exit(1)
+print("All stats tests passed")
diff --git a/tests/test_stats_report.py b/tests/test_stats_report.py
new file mode 100755
index 0000000..9ca9f30
--- /dev/null
+++ b/tests/test_stats_report.py
@@ -0,0 +1,228 @@
+#!/usr/bin/env python3
+# SPDX-License-Identifier: MIT
+#
+# Checks dto-stats-report's arithmetic on made-up stats runs with known
+# answers, its perf sample categories on made-up "perf report" lines, and
+# that it finds libc's mem* functions by address in this machine's libc.
+# Needs neither DSA nor perf.
+#
+# Usage: test_stats_report.py <dto-stats-report>
+
+import ctypes
+import importlib.machinery
+import importlib.util
+import json
+import os
+import subprocess
+import sys
+import tempfile
+
+REPORT = sys.argv[1]
+failures = 0
+sys.dont_write_bytecode = True      # importing the report must not litter
+
+
+def check(desc, ok, detail=""):
+    global failures
+    print(("PASS: " if ok else "FAIL: ") + desc + ("" if ok else "  " + str(detail)))
+    failures += not ok
+
+
+def close(a, b):
+    return a is not None and abs(a - b) <= 1e-6 * max(1, abs(b))
+
+
+def ops(**kw):
+    """ops(memcpy_cpu=(calls, bytes, ns), ...) in the summary's layout"""
+    out = {}
+    for key, (calls, nbytes, ns) in kw.items():
+        op, path = key.split("_", 1)
+        out.setdefault(op, {})[path] = {"calls": calls, "bytes": nbytes, "ns": ns}
+    return out
+
+
+def summary(pid, totals, hist, stdc_only=0, page_faults=0):
+    return {"type": "summary", "pid": pid, "program": "app", "run_ms": 1.0,
+            "config": {"stdc_only": stdc_only, "wait_method": "busypoll",
+                       "min_bytes": 65536, "cpu_size_fraction": 0.0,
+                       "auto_adjust": 1, "numa_aware": 0, "wqs": 1},
+            "totals": totals,
+            "failures": {"retry": 0, "page_fault": page_faults, "other": 0},
+            "histogram": hist}
+
+
+def write_run_dir(runs):
+    """runs: (name, mode, wall_s, cpu_s, [summary records], exit_status)"""
+    d = tempfile.mkdtemp()
+    with open(os.path.join(d, "runs.jsonl"), "w") as f:
+        f.write(json.dumps({"type": "stats_run", "library": "/x/libdto.so",
+                            "host": "test", "repeat": 1, "perf": None,
+                            "command": ["app"], "env": {}}) + "\n")
+        for name, mode, wall, cpu, recs, status in runs:
+            f.write(json.dumps({"type": "run", "name": name, "mode": mode,
+                                "round": 1, "profiled": False,
+                                "stats": name + ".jsonl", "wall_s": wall,
+                                "user_s": cpu, "sys_s": 0.0, "maxrss_kb": 1,
+                                "exit_status": status}) + "\n")
+            with open(os.path.join(d, name + ".jsonl"), "w") as s:
+                for r in recs:
+                    s.write(json.dumps(r) + "\n")
+    return d
+
+
+def report(d):
+    p = subprocess.run([REPORT, "--json", d], capture_output=True, text=True)
+    if p.returncode:
+        print(p.stderr)
+        return None
+    return json.loads(p.stdout)
+
+
+MiB = 1 << 20
+B_CPY, B_SET = (1024 * MiB, 512 * MiB)
+H_CPY = (1 * MiB, 1 * MiB + MiB // 4 - 1)
+H_SET = (256 << 10, (320 << 10) - 1)
+
+
+def hist(cpy, mset):
+    return [{"min": H_SET[0], "max": H_SET[1], "ops": mset},
+            {"min": H_CPY[0], "max": H_CPY[1], "ops": cpy}]
+
+
+# CPU-only: memcpy 4 s and memset 1 s of 10 s CPU time, so f = 0.5
+cpu_cpy = ops(memcpy_cpu=(1000, B_CPY, 4e9))
+cpu_set = ops(memset_cpu=(1000, B_SET, 1e9))
+cpu_rec = summary(1, {**cpu_cpy, **cpu_set}, hist(cpu_cpy, cpu_set))
+# the same work in a run whose threads were preempted inside the calls:
+# three times the mem* time, a little more CPU time
+slow_cpy = ops(memcpy_cpu=(1000, B_CPY, 12e9))
+slow_set = ops(memset_cpu=(1000, B_SET, 3e9))
+slow_rec = summary(1, {**slow_cpy, **slow_set}, hist(slow_cpy, slow_set))
+# with DSA: memcpy 4x faster; memset 90% of bytes on DSA, 0.3 s in all.
+# Two processes, each reporting part of the work.
+dsa_cpy = ops(memcpy_dsa=(1000, B_CPY, 1e9))
+dsa_set = ops(memset_dsa=(900, B_SET * 9 // 10, 2e8),
+              memset_cpu=(100, B_SET // 10, 1e8))
+dsa_recs = [summary(2, dsa_cpy, hist(dsa_cpy, {}), page_faults=5),
+            summary(3, dsa_set, hist({}, dsa_set))]
+
+# --- the numbers ---
+d = write_run_dir([
+    ("cpu-1", "cpu", 10.0, 10.0, [cpu_rec], 0),
+    ("dsa-1", "dsa", 6.0, 6.0, dsa_recs, 0),
+    ("cpu-2", "cpu", 30.0, 10.5, [slow_rec], 0),    # an outlier
+    ("dsa-2", "dsa", 7.0, 6.2, dsa_recs, 0),
+    ("cpu-3", "cpu", 11.0, 10.0, [cpu_rec], 0),
+    ("dsa-3", "dsa", 6.5, 6.1, dsa_recs, 0),
+])
+r = report(d)
+check("report runs", r is not None)
+if r:
+    a = r["amdahl"]
+    s = 5.0 / 1.3
+    check("f = mem* share of CPU-only CPU time", close(a["f"], 0.5), a["f"])
+    check("s = mem* speedup", close(a["s"], s), a["s"])
+    check("predicted = 1 / ((1 - f) + f / s)", close(a["predicted"], 1 / (0.5 + 0.5 / s)),
+          a["predicted"])
+    check("limit = 1 / (1 - f)", close(a["limit"], 2.0), a["limit"])
+    check("measured speedup uses median wall times", close(r["measured"]["wall"], 11 / 6.5),
+          r["measured"])
+    check("median CPU time", close(r["modes"]["dsa"]["cpu_s"], 6.1), r["modes"]["dsa"])
+    check("mem* time is the median-CPU run's, not a mean over runs",
+          close(r["modes"]["cpu"]["mem_s"], 5.0), r["modes"]["cpu"]["mem_s"])
+    check("time outside mem* calls comes from that same run",
+          close(r["outside_mem_s"]["cpu"], 5.0), r["outside_mem_s"])
+    check("processes of a run are summed", close(r["modes"]["dsa"]["mem_s"], 1.3),
+          r["modes"]["dsa"]["mem_s"])
+    o = r["ops"]
+    check("memcpy speedup", close(o["memcpy"]["speedup"], 4.0), o["memcpy"])
+    check("memset speedup", close(o["memset"]["speedup"], 1 / 0.3), o["memset"])
+    check("memset share on DSA", close(o["memset"]["on_dsa"], 0.9), o["memset"])
+    check("memcpy share of CPU time", close(o["memcpy"]["share"], 0.4), o["memcpy"])
+    check("memset calls", o["memset"]["calls"] == 1000, o["memset"])
+    check("page faults per run", close(r["dsa_failures"]["page_fault"], 5),
+          r["dsa_failures"])
+    b = {x["min"]: x for x in r["buckets"]}
+    check("memcpy size range: average call time",
+          close(b[H_CPY[0]]["cpu_only_us"], 4000) and close(b[H_CPY[0]]["with_dsa_us"], 1000),
+          b.get(H_CPY[0]))
+    check("no notes for a clean run", r["notes"] == [], r["notes"])
+
+# --- the DSA run used no DSA ---
+fallback = summary(2, {**cpu_cpy, **cpu_set}, hist(cpu_cpy, cpu_set), stdc_only=1)
+r = report(write_run_dir([("cpu-1", "cpu", 10.0, 10.0, [cpu_rec], 0),
+                          ("dsa-1", "dsa", 10.0, 10.0, [fallback], 0)]))
+check("says when DTO found no work queue",
+      r and any("no usable work queue" in n for n in r["notes"]), r and r["notes"])
+
+# --- different amounts of work ---
+more = ops(memcpy_dsa=(2000, 2 * B_CPY, 2e9))
+r = report(write_run_dir([("cpu-1", "cpu", 10.0, 10.0, [cpu_rec], 0),
+                          ("dsa-1", "dsa", 10.0, 10.0,
+                           [summary(2, {**more, **dsa_set}, [])], 0)]))
+check("says when the runs did different work",
+      r and any("different amounts" in n for n in r["notes"]), r and r["notes"])
+check("speedup is per byte", r and close(r["ops"]["memcpy"]["speedup"], 4.0),
+      r and r["ops"]["memcpy"])
+
+# --- DSA slower for small calls, faster for big ones ---
+small_c, small_d = ops(memset_cpu=(1000, B_SET, 1e8)), ops(memset_dsa=(1000, B_SET, 3e8))
+r = report(write_run_dir([
+    ("cpu-1", "cpu", 10.0, 10.0, [summary(1, {}, hist(cpu_cpy, small_c))], 0),
+    ("dsa-1", "dsa", 10.0, 10.0, [summary(2, {}, hist(dsa_cpy, small_d))], 0)]))
+check("suggests DTO_MIN_BYTES from where DSA starts winning",
+      r and any(f"DTO_MIN_BYTES={H_CPY[0]}" in n for n in r["notes"]), r and r["notes"])
+
+# --- a failed run ---
+r = report(write_run_dir([("cpu-1", "cpu", 10.0, 10.0, [cpu_rec], 0),
+                          ("dsa-1", "dsa", 1.0, 1.0, [], 3)]))
+check("says a run failed", r and any("exited with status 3" in n for n in r["notes"]),
+      r and r["notes"])
+
+p = subprocess.run([REPORT, tempfile.mkdtemp()], capture_output=True, text=True)
+check("a directory without runs.jsonl is an error", p.returncode != 0, p.stdout)
+
+# --- perf sample categories ---
+loader = importlib.machinery.SourceFileLoader("dto_stats_report", REPORT)
+spec = importlib.util.spec_from_loader("dto_stats_report", loader)
+rep = importlib.util.module_from_spec(spec)
+loader.exec_module(rep)
+lines = ["400|libc.so.6           |[.] __memmove_avx512_unaligned_erms",
+         "100|libc.so.6           |[.] __memset_evex_unaligned_erms",
+         "100|libc.so.6           |[.] __memcmp_evex_movbe",
+         "100|app                 |[.] crc32c_hw",
+         "100|[kernel.kallsyms]   |[k] rep_movs_alternative",
+         "50|[kernel.kallsyms]    |[k] clear_page_erms",
+         "50|[kernel.kallsyms]    |[k] schedule",
+         "50|libdto.so            |[.] dsa_wait_busy_poll",
+         "40|app                  |[.] main",
+         "10|libc.so.6            |[.] 0x00000000000fffff",
+         "not a sample line"]
+shares, unknown = rep.classify(lines)
+want = {"memcpy/memmove": 0.4, "memset": 0.1, "memcmp": 0.1, "CRC": 0.1,
+        "kernel copy/clear": 0.15, "other kernel code": 0.05, "DTO library": 0.05,
+        "other user code": 0.04, "libc, unknown function": 0.01}
+check("perf samples by category", all(close(shares[k], v) for k, v in want.items()), shares)
+check("unplaced libc share", close(unknown, 0.01), unknown)
+
+# --- libc functions by address, in this machine's libc ---
+loc = rep.LibcLocator()
+if not loc.ok:
+    print("SKIP: libc located by address (no build ID or not 64-bit ELF)")
+else:
+    libc = ctypes.CDLL("libc.so.6")
+    base = [int(line.split("-")[0], 16) for line in open("/proc/self/maps")
+            if rep.is_libc(line.split()[-1]) and int(line.split()[2], 16) == 0][0]
+    for name, cat in (("memmove", "memcpy/memmove"), ("memset", "memset"),
+                      ("memcmp", "memcmp")):
+        vaddr = ctypes.cast(getattr(libc, name), ctypes.c_void_p).value - base + 16
+        # perf reports file offsets
+        off = [vaddr - v + o for o, v, n in loc.loads if v <= vaddr < v + n][0]
+        check(f"a sample just inside {name} is {cat}", loc.category(off) == cat,
+              loc.category(off))
+
+print()
+if failures:
+    print(f"{failures} stats report test(s) failed")
+    sys.exit(1)
+print("All stats report tests passed")
diff --git a/tests/test_stats_workload.c b/tests/test_stats_workload.c
new file mode 100644
index 0000000..d6f257a
--- /dev/null
+++ b/tests/test_stats_workload.c
@@ -0,0 +1,132 @@
+/*******************************************************************************
+ * Copyright (C) 2026 Intel Corporation
+ *
+ * SPDX-License-Identifier: MIT
+ *
+ * test_stats_workload.c - known memory operations for tests/test_stats.py,
+ * which runs this with DTO_COLLECT_STATS set and checks what DTO reports.
+ *
+ *   basic <threads> <iters> [sleep_ms]  each thread, per iteration: memcpy,
+ *                                       memset, memmove, memcmp of
+ *                                       SIZE_{CPY,SET,MOV,CMP} bytes
+ *   churn <threads> <iters>             the same, from short-lived threads
+ *                                       started one after another; prints
+ *                                       the mapping count before and after
+ *   fork <iters>                        memcpy SIZE_CPY in the parent and
+ *                                       SIZE_CHILD in a forked child
+ *
+ * Built with -fno-builtin so every call reaches DTO.
+ ******************************************************************************/
+
+#include <pthread.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <sys/mman.h>
+#include <sys/wait.h>
+#include <time.h>
+#include <unistd.h>
+
+/* distinctive sizes, each in its own histogram bucket */
+#define SIZE_CPY (131072 + 123)
+#define SIZE_SET (65536 + 4096)
+#define SIZE_MOV (98304 + 7)
+#define SIZE_CMP (262144 + 64)
+#define SIZE_CHILD (524288 + 1)
+#define BUF (1 << 20)
+
+static int iters, sleep_ms;
+
+static void *work(void *arg)
+{
+	char *a = mmap(NULL, BUF, PROT_READ | PROT_WRITE,
+		       MAP_PRIVATE | MAP_ANONYMOUS | MAP_POPULATE, -1, 0);
+	char *b = mmap(NULL, BUF, PROT_READ | PROT_WRITE,
+		       MAP_PRIVATE | MAP_ANONYMOUS | MAP_POPULATE, -1, 0);
+	struct timespec pause = { 0, sleep_ms * 1000000L };
+
+	for (int i = 0; i < iters; i++) {
+		memcpy(b, a, SIZE_CPY);
+		memset(a, i, SIZE_SET);
+		memmove(b + 64, b, SIZE_MOV);
+		if (memcmp(a, a + 0, SIZE_CMP))
+			abort();
+		if (sleep_ms)
+			nanosleep(&pause, NULL);
+	}
+	munmap(a, BUF);
+	munmap(b, BUF);
+	return NULL;
+}
+
+static int count_maps(void)
+{
+	char line[512];
+	int n = 0;
+	FILE *f = fopen("/proc/self/maps", "r");
+
+	while (f && fgets(line, sizeof(line), f))
+		n++;
+	if (f)
+		fclose(f);
+	return n;
+}
+
+static void copies(size_t size, int n)
+{
+	char *a = mmap(NULL, BUF, PROT_READ | PROT_WRITE,
+		       MAP_PRIVATE | MAP_ANONYMOUS | MAP_POPULATE, -1, 0);
+	char *b = mmap(NULL, BUF, PROT_READ | PROT_WRITE,
+		       MAP_PRIVATE | MAP_ANONYMOUS | MAP_POPULATE, -1, 0);
+
+	for (int i = 0; i < n; i++)
+		memcpy(b, a, size);
+	munmap(a, BUF);
+	munmap(b, BUF);
+}
+
+int main(int argc, char **argv)
+{
+	pthread_t t[256];
+
+	if (argc < 3)
+		return 2;
+	if (!strcmp(argv[1], "basic")) {
+		int nt = atoi(argv[2]);
+
+		iters = atoi(argv[3]);
+		sleep_ms = argc > 4 ? atoi(argv[4]) : 0;
+		for (int i = 0; i < nt; i++)
+			pthread_create(&t[i], NULL, work, NULL);
+		for (int i = 0; i < nt; i++)
+			pthread_join(t[i], NULL);
+	} else if (!strcmp(argv[1], "churn")) {
+		int nt = atoi(argv[2]), before, after;
+
+		iters = atoi(argv[3]);
+		/* the first thread's block and allocations, so they don't count */
+		pthread_create(&t[0], NULL, work, NULL);
+		pthread_join(t[0], NULL);
+		before = count_maps();
+		for (int i = 1; i < nt; i++) {
+			pthread_create(&t[0], NULL, work, NULL);
+			pthread_join(t[0], NULL);
+		}
+		after = count_maps();
+		printf("maps_before %d maps_after %d\n", before, after);
+	} else if (!strcmp(argv[1], "fork")) {
+		int n = atoi(argv[2]), st;
+		pid_t pid;
+
+		copies(SIZE_CPY, n);
+		pid = fork();
+		if (pid == 0) {
+			copies(SIZE_CHILD, n);
+			exit(0);
+		}
+		waitpid(pid, &st, 0);
+		copies(SIZE_CPY, n);
+		printf("child %d\n", pid);
+	}
+	return 0;
+}