diff --git a/README.md b/README.md
index b8f2b58..2e3701b 100644
--- a/README.md
+++ b/README.md
@@ -163,9 +163,26 @@ skim_dir: /ceph/USER/shapesmith/TAG
 output_dir: output/TAG
 ml_dir: /ceph/USER/shapesmith_ml/TAG
 workers: 16
+log_level: INFO                           # DEBUG | INFO | WARNING | ERROR
 combine: {cmssw_dir: /work/USER/CMSSW_14_1_9, scram_arch: el9_amd64_gcc12}
 ```
 
+## Logging
+
+`log_level` sets the level of ShapeSmith and of the analysis package; one run changes it with
+`-s log_level=DEBUG`. Other libraries (matplotlib, numexpr, uproot, ...) log from WARNING on.
+
+- Console: the messages, on stderr (stdout keeps the reports of `validate` and `inspect`).
+- Log file: every command but `validate` and `inspect` also writes
+  `<output_dir>/logs/<command>_<YYYYmmdd-HHMMSS>.log` with time, level, process and logger, and the
+  traceback of a failing command.
+- INFO: the run (versions, configuration, channels, workers), the analysis per channel, what every
+  stage does per channel, the progress of the parallel jobs (every tenth), the files written and
+  the duration. DEBUG adds every ntuple file, booking, estimate, datacard bin, fit and plot, and
+  the output of the CMSSW commands (at INFO it is shown only when they fail).
+- The worker processes log through the main process, so their messages reach the console and
+  the log file as well.
+
 ## Notes
 
 - Expressions (cuts, weights, variables) are pandas `eval` syntax on ntuple columns, evaluated
diff --git a/src/shapesmith/binning.py b/src/shapesmith/binning.py
index b44b155..c5cd7f4 100644
--- a/src/shapesmith/binning.py
+++ b/src/shapesmith/binning.py
@@ -12,6 +12,7 @@
 
 import dataclasses
 import json
+import logging
 from pathlib import Path
 from typing import TYPE_CHECKING
 
@@ -28,6 +29,8 @@
 EPSILON = 1e-4
 RECORD = "binning.json"
 
+logger = logging.getLogger(__name__)
+
 
 def inside(values: np.ndarray, rule: EqualData) -> np.ndarray:
     """The finite values in (low, high), as float64."""
@@ -68,7 +71,9 @@ def resolve(config: RunConfig, analysis: Analysis, channel_name: str, targets: l
         except ValueError as error:
             raise ValueError(f"{channel_name}/{target.category}/{target.variable.name}: equal-data binning: {error}") from None
         resolved.append(dataclasses.replace(target, variable=dataclasses.replace(target.variable, edges=edges)))
-        record.setdefault(target.category, {})[target.variable.name] = {"rule": dataclasses.asdict(rule), "n_data": int(inside(values, rule).size), "edges": list(edges)}
+        n_data = int(inside(values, rule).size)
+        record.setdefault(target.category, {})[target.variable.name] = {"rule": dataclasses.asdict(rule), "n_data": n_data, "edges": list(edges)}
+        logger.info(f"{channel_name}/{target.category}/{target.variable.name}: {len(edges) - 1} equal-data bins from {n_data} data events, edges {', '.join(f'{edge:.6g}' for edge in edges)}")
     return resolved, record
 
 
diff --git a/src/shapesmith/cli.py b/src/shapesmith/cli.py
index 9301b89..5dc8557 100644
--- a/src/shapesmith/cli.py
+++ b/src/shapesmith/cli.py
@@ -1,37 +1,85 @@
 """Command-line interface: one sub-command per analysis step."""
 from __future__ import annotations
 
+import functools
+import json
 import logging
+import time
 from pathlib import Path
 from typing import Optional
 
 import typer
 from rich.console import Console
-from rich.logging import RichHandler
 from rich.table import Table
 
-from shapesmith import __version__
+from shapesmith import __version__, logs
 from shapesmith.config import RunConfig, load_analysis, load_config
-from shapesmith.provenance import write_versions
+from shapesmith.ntuples import join_url
+from shapesmith.provenance import record, write_versions
 
 app = typer.Typer(help=f"ShapeSmith {__version__}: CROWN ntuples -> skims -> histograms -> datacards, fits, plots, ML folds, measurements.", no_args_is_help=True)
+logger = logging.getLogger("shapesmith.cli")  # also under `python -m shapesmith.cli`, where __name__ is __main__
+_command = {"name": "shapesmith"}  # the running sub-command, set by _logged (log file name and first log line)
 
 ConfigOption = typer.Option(..., "--config", "-c", help="run configuration YAML")
 ChannelsOption = typer.Option(None, "--channels", help="comma separated subset of the configured channels")
 WorkersOption = typer.Option(None, "--workers", help="override `workers` of the configuration (1 = inline, best for debugging)")
-SetOption = typer.Option(None, "--set", "-s", help="override a configuration entry, e.g. -s switches.jet_fakes=ff -s workers=4 (repeatable, values parsed as YAML)")
-
-
-def _setup(config_path: Path, channels: Optional[str], log_level: str = "INFO", workers: Optional[int] = None, overrides: Optional[list[str]] = None):
-    logging.basicConfig(level=log_level, format="%(message)s", handlers=[RichHandler(show_path=False)], force=True)
+SetOption = typer.Option(None, "--set", "-s", help="override a configuration entry, e.g. -s switches.jet_fakes=ff -s workers=4 -s log_level=DEBUG (repeatable, values parsed as YAML)")
+
+
+def _logged(command):
+    """Log the end of a command: its duration, or its error (with the traceback in the log file); then close the log."""
+
+    @functools.wraps(command)
+    def wrapper(*args, **kwargs):
+        _command["name"] = command.__name__.replace("_", "-")
+        start = time.monotonic()
+        try:
+            result = command(*args, **kwargs)
+        except Exception as error:
+            if logs.active():
+                logger.error(f"failed after {logs.duration(time.monotonic() - start)}: {type(error).__name__}: {error}", exc_info=True)
+            raise
+        else:
+            if logs.active():
+                logger.info(f"done in {logs.duration(time.monotonic() - start)}")
+            return result
+        finally:
+            logs.reset()
+
+    return wrapper
+
+
+def _git(entry: dict | None) -> str:
+    return f" (git {entry['git'][:10]})" if entry and entry.get("git") else ""
+
+
+def _log_start(config: RunConfig, config_path: Path, overrides: list[str], selected: list[str], log_file: Path | None) -> None:
+    """What runs: command, configuration, overrides, analysis, versions and git hashes (provenance.record)."""
+    versions = record(config)
+    logger.info(f"ShapeSmith {__version__} {_command['name']}: {Path(config_path).absolute()}" + (f", overrides {' '.join(overrides)}" if overrides else ""))
+    logger.info(f"analysis {config.analysis}{_git(versions['analysis'])}, era {config.era}, channels {','.join(selected)}, {config.workers} workers, log level {config.log_level}")
+    database = f", sample database {config.sample_database}{_git(versions['sample_database'])}" if config.sample_database else ""
+    logger.info(f"shapesmith{_git(versions['shapesmith'])}{database}, ntuples {join_url(config.ntuples.server, config.ntuples.base)}")
+    if log_file is not None:
+        logger.info(f"log file {log_file}")
+    logger.debug(f"resolved configuration: {json.dumps(versions['config'], sort_keys=True)}")
+
+
+def _setup(config_path: Path, channels: Optional[str], workers: Optional[int] = None, overrides: Optional[list[str]] = None, log_file: bool = True):
+    """The configuration (+ overrides), the logging of the command (console, and <output_dir>/logs/ unless `log_file`
+    is False) and the validated analysis with the selected channels."""
     config = load_config(config_path, overrides or [])
     if workers is not None:
         config = config.model_copy(update={"workers": workers})
-    analysis = load_analysis(config)
+    path = logs.log_path(config.output_dir, _command["name"]) if log_file else None
+    logs.configure(config.log_level, path, [config.analysis.partition(":")[0].split(".")[0]])
     selected = channels.split(",") if channels else list(config.channels)
+    _log_start(config, config_path, overrides or [], selected, path)
     unknown = set(selected) - set(config.channels)
     if unknown:
         raise typer.BadParameter(f"channels {sorted(unknown)} are not in the configuration")
+    analysis = load_analysis(config)
     return config, analysis, selected
 
 
@@ -51,11 +99,12 @@ def _split(value: Optional[str]) -> Optional[list[str]]:
 
 
 @app.command()
+@_logged
 def validate(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption):
     """Load and validate the analysis; report processes, selections and the columns every channel needs."""
     from shapesmith.skim import needed_columns
 
-    cfg, analysis, selected = _setup(config, channels, overrides=overrides)
+    cfg, analysis, selected = _setup(config, channels, overrides=overrides, log_file=False)
     console = Console()
     console.print(f"[bold]{analysis.name}[/bold]  era {analysis.era}, {analysis.lumi_pb:g} pb^-1, switches {dict(cfg.switches)}")
     for name in selected:
@@ -77,6 +126,7 @@ def validate(config: Path = ConfigOption, channels: Optional[str] = ChannelsOpti
 
 
 @app.command()
+@_logged
 def skim(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, force: bool = typer.Option(False, "--force", help="skim again, ignoring stored skims"), samples: Optional[str] = typer.Option(None, "--samples", help="comma separated sample groups or nick prefixes (default: all)"), workers: Optional[int] = WorkersOption):
     """Stage 1: ntuples (+ friends) -> Parquet skims."""
     from shapesmith.skim import run_skim
@@ -84,10 +134,11 @@ def skim(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption,
     cfg, analysis, selected = _setup(config, channels, workers=workers, overrides=overrides)
     write_versions(cfg, analysis.name, cfg.skim_dir)
     results = run_skim(cfg, analysis, selected, force=force, samples=_split(samples))
-    typer.echo(f"{sum(1 for r in results if not r.skipped)} files skimmed, {sum(1 for r in results if r.skipped)} reused, {sum(r.n_out for r in results)} events kept")
+    logger.info(f"{sum(1 for r in results if not r.skipped)} files skimmed, {sum(1 for r in results if r.skipped)} reused, {sum(r.n_out for r in results)} events kept")
 
 
 @app.command()
+@_logged
 def hist(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, control: bool = typer.Option(False, "--control", help="control variables instead of NN categories"), variables: Optional[str] = typer.Option(None, help="comma separated control variables"), regions: Optional[str] = typer.Option(None, "--regions", help="comma separated regions, or 'all' (default: nominal plus estimator regions)"), skip_systematics: bool = typer.Option(False, "--skip-systematics", help="no weight or column variations"), processes: Optional[str] = typer.Option(None, help="comma separated process names"), workers: Optional[int] = WorkersOption):
     """Stage 2: skims -> histograms (control_shapes.root or shapes.root); other histograms in the file are kept."""
     from shapesmith.fill import run_hist
@@ -95,11 +146,11 @@ def hist(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption,
     cfg, analysis, selected = _setup(config, channels, workers=workers, overrides=overrides)
     write_versions(cfg, analysis.name, cfg.output_dir)
     output = _shapes_path(cfg, control)
-    hset = run_hist(cfg, analysis, selected, control, _split(variables), not skip_systematics, _split(processes), output, _split(regions))
-    typer.echo(f"{len(hset)} histograms in {output}")
+    run_hist(cfg, analysis, selected, control, _split(variables), not skip_systematics, _split(processes), output, _split(regions))
 
 
 @app.command()
+@_logged
 def estimate(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, control: bool = typer.Option(False, "--control")):
     """Stage 3: add the estimated processes and variations of every channel to the histogram file."""
     from shapesmith.estimates import run_estimates
@@ -110,19 +161,21 @@ def estimate(config: Path = ConfigOption, channels: Optional[str] = ChannelsOpti
     hset = HistogramSet.load(path)
     added = run_estimates(hset, analysis, selected)
     hset.save(path)
-    typer.echo(f"{len(added)} histograms added to {path}")
+    logger.info(f"{len(added)} histograms added to {path}")
 
 
 @app.command()
-def measure(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, suggest_binning: bool = typer.Option(False, "--suggest-binning", help="print proposed bin edges instead of measuring"), merge: bool = typer.Option(False, "--merge", help="combine the results of the earlier runs (e.g. one per working point) into one payload"), workers: Optional[int] = WorkersOption):
+@_logged
+def measure(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, suggest_binning: bool = typer.Option(False, "--suggest-binning", help="log proposed bin edges instead of measuring"), merge: bool = typer.Option(False, "--merge", help="combine the results of the earlier runs (e.g. one per working point) into one payload"), workers: Optional[int] = WorkersOption):
     """Run the analysis's measurement; results in <output_dir>/<measurement>/<era>/."""
     from shapesmith.measurements import run_measure
 
     cfg, analysis, selected = _setup(config, channels, workers=workers, overrides=overrides)
-    typer.echo(f"measurement output: {run_measure(cfg, analysis, selected, suggest_binning, merge)}")
+    logger.info(f"measurement output: {run_measure(cfg, analysis, selected, suggest_binning, merge)}")
 
 
 @app.command()
+@_logged
 def sync(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption):
     """Write combine-style shape files per channel."""
     from shapesmith.histogram import HistogramSet
@@ -130,10 +183,11 @@ def sync(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption,
 
     cfg, analysis, selected = _setup(config, channels, overrides=overrides)
     for path in run_sync(HistogramSet.load(_shapes_path(cfg, False)), analysis, selected, cfg.output_dir):
-        typer.echo(f"wrote {path}")
+        logger.info(f"wrote {path}")
 
 
 @app.command()
+@_logged
 def datacards(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, systematics: bool = typer.Option(True, "--systematics/--no-systematics"), min_background: Optional[float] = typer.Option(1.0, help="merge bins below this background yield; negative disables")):
     """Write text datacards and their shapes file per final state (each channel and the combination)."""
     from shapesmith.datacards import run_datacards
@@ -142,10 +196,11 @@ def datacards(config: Path = ConfigOption, channels: Optional[str] = ChannelsOpt
     cfg, analysis, selected = _setup(config, channels, overrides=overrides)
     threshold = None if min_background is None or min_background < 0 else min_background
     for name, path in run_datacards(HistogramSet.load(_shapes_path(cfg, False)), analysis, _final_states(selected), cfg.output_dir / "datacards", systematics, threshold).items():
-        typer.echo(f"{name}: {path}")
+        logger.info(f"datacard {name}: {path}")
 
 
 @app.command()
+@_logged
 def fit(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, final_states: Optional[str] = typer.Option(None, help="comma separated; default: each channel and 'all'"), skip_combine: bool = typer.Option(False, "--skip-combine", help="only collect existing combine outputs")):
     """Run combine (expected limits, significance, best fit) and summarise the results."""
     from shapesmith.limits import run_limits
@@ -154,10 +209,11 @@ def fit(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, o
     states = _split(final_states) or list(_final_states(selected))
     results = run_limits(cfg, cfg.output_dir / "datacards", states, skip_combine)
     typer.echo((cfg.output_dir / "datacards" / "limits.md").read_text())
-    typer.echo(f"{len(results)} final states summarised")
+    logger.info(f"{len(results)} final states summarised")
 
 
 @app.command()
+@_logged
 def plot(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, control: bool = typer.Option(False, "--control"), category: Optional[str] = typer.Option(None), variables: Optional[str] = typer.Option(None), region: str = typer.Option("nominal", "--region", help="histogram region to plot"), blind: bool = typer.Option(False, "--blind"), log: bool = typer.Option(False, "--log"), signal_scale: Optional[float] = typer.Option(None), normalize_by_bin_width: bool = typer.Option(False, "--normalize-by-bin-width")):
     """Prefit stack plots of control variables (--control) or NN categories."""
     from shapesmith.histogram import HistogramSet
@@ -165,10 +221,11 @@ def plot(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption,
 
     cfg, analysis, selected = _setup(config, channels, overrides=overrides)
     files = run_plot(HistogramSet.load(_shapes_path(cfg, control)), analysis, selected, control, category, _split(variables), cfg.output_dir / "plots", region=region, blind=blind, log=log, signal_scale=signal_scale, normalize_by_bin_width=normalize_by_bin_width)
-    typer.echo(f"{len(files)} files written to {cfg.output_dir / 'plots'}")
+    logger.info(f"{len(files)} files written to {cfg.output_dir / 'plots'}")
 
 
 @app.command("ml-export")
+@_logged
 def ml_export(config: Path = ConfigOption, channels: Optional[str] = ChannelsOption, overrides: Optional[list[str]] = SetOption, workers: Optional[int] = WorkersOption):
     """Training folds (Feather) from the skims."""
     from shapesmith.ml_export import run_ml_export
@@ -177,7 +234,7 @@ def ml_export(config: Path = ConfigOption, channels: Optional[str] = ChannelsOpt
     if cfg.ml_dir is not None:
         write_versions(cfg, analysis.name, cfg.ml_dir)
     for path in run_ml_export(cfg, analysis, selected):
-        typer.echo(f"wrote {path}")
+        logger.info(f"wrote {path}")
 
 
 @app.command()
diff --git a/src/shapesmith/cmssw.py b/src/shapesmith/cmssw.py
index ca0d0d6..d44d6ed 100644
--- a/src/shapesmith/cmssw.py
+++ b/src/shapesmith/cmssw.py
@@ -2,7 +2,8 @@
 
 The subshell starts from a clean environment: the Python environment ShapeSmith runs in (e.g. an LCG view) would
 otherwise shadow CMSSW's libraries and Python packages. With a timeout, the whole process group of the subshell is
-killed when it expires (a fit can iterate without end).
+killed when it expires (a fit can iterate without end). Its output is logged line by line at DEBUG, and its last
+lines at ERROR when it fails.
 """
 from __future__ import annotations
 
@@ -10,6 +11,8 @@
 import os
 import signal
 import subprocess
+import threading
+from collections import deque
 from pathlib import Path
 from typing import Sequence
 
@@ -17,6 +20,8 @@
 
 logger = logging.getLogger(__name__)
 KEPT_VARIABLES = ("HOME", "USER", "LOGNAME", "TMPDIR", "X509_USER_PROXY")  # the environment a CMSSW subshell starts with, besides PATH
+TAIL = 40  # lines of output logged with a failure
+READER_TIMEOUT = 10.0  # seconds to wait for the rest of the output after the subshell ended
 
 
 def script(combine: CombineConfig, cwd: Path, commands: Sequence[str]) -> str:
@@ -38,18 +43,39 @@ def clean_environment() -> dict[str, str]:
     return {"PATH": "/usr/bin:/bin", **{name: os.environ[name] for name in KEPT_VARIABLES if name in os.environ}}
 
 
+def _read(stream, tail: deque) -> None:
+    for line in stream:
+        tail.append(line.rstrip("\n"))
+        logger.debug(f"cmssw: {tail[-1]}")
+    stream.close()
+
+
+def _last_lines(tail: deque) -> str:
+    return ", last output:\n" + "\n".join(tail) if tail else ", no output"
+
+
 def run(commands: Sequence[str], combine: CombineConfig | None, cwd: Path, timeout: float | None = None) -> None:
     """Raise CalledProcessError on a failing command, TimeoutExpired after `timeout` seconds."""
     if combine is None:
         raise ValueError("the run configuration needs `combine` (cmssw_dir) to run CMSSW commands")
     logger.info(f"running {len(commands)} CMSSW commands in {cwd}")
-    arguments = ["bash", "-c", script(combine, cwd, commands)]
-    process = subprocess.Popen(arguments, env=clean_environment(), start_new_session=True)
+    text = script(combine, cwd, commands)
+    logger.debug(f"CMSSW script:\n{text}")
+    arguments = ["bash", "-c", text]
+    process = subprocess.Popen(arguments, env=clean_environment(), start_new_session=True, stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True, errors="replace")
+    tail: deque = deque(maxlen=TAIL)
+    reader = threading.Thread(target=_read, args=(process.stdout, tail), daemon=True)
+    reader.start()
     try:
         code = process.wait(timeout=timeout)
     except subprocess.TimeoutExpired:
         os.killpg(process.pid, signal.SIGKILL)
         process.wait()
+        reader.join(READER_TIMEOUT)
+        logger.error(f"CMSSW commands in {cwd} killed after {timeout:g} s{_last_lines(tail)}")
         raise
+    reader.join(READER_TIMEOUT)
     if code:
+        logger.error(f"CMSSW commands in {cwd} failed with exit code {code}{_last_lines(tail)}")
         raise subprocess.CalledProcessError(code, arguments)
+    logger.debug(f"CMSSW commands in {cwd} done")
diff --git a/src/shapesmith/config.py b/src/shapesmith/config.py
index 7182d53..64a2094 100644
--- a/src/shapesmith/config.py
+++ b/src/shapesmith/config.py
@@ -2,16 +2,20 @@
 from __future__ import annotations
 
 import importlib
+import logging
 import os
+from collections import Counter
 from pathlib import Path
-from typing import Any, Iterable
+from typing import Any, Iterable, Literal
 
 import yaml
-from pydantic import BaseModel, ConfigDict, Field
+from pydantic import BaseModel, ConfigDict, Field, field_validator
 
 from shapesmith.model import Analysis, SampleKind
 from shapesmith.validate import validate
 
+logger = logging.getLogger(__name__)
+
 
 class FriendConfig(BaseModel):
     model_config = ConfigDict(extra="forbid")
@@ -44,8 +48,14 @@ class RunConfig(BaseModel):
     ml_dir: Path | None = None
     sample_database: Path | None = None  # KingMaker datasets.json used for the normalisation
     workers: int = 4
+    log_level: Literal["DEBUG", "INFO", "WARNING", "ERROR"] = "INFO"  # of ShapeSmith and the analysis (logs.configure)
     combine: CombineConfig | None = None
 
+    @field_validator("log_level", mode="before")
+    @classmethod
+    def _upper_case(cls, value: Any) -> Any:
+        return value.upper() if isinstance(value, str) else value
+
 
 RELATIVE_PATH_FIELDS = ("skim_dir", "output_dir", "ml_dir", "sample_database")
 
@@ -77,6 +87,20 @@ def load_config(path: str | Path, overrides: Iterable[str] = ()) -> RunConfig:
     return RunConfig.model_validate(raw)
 
 
+def summary(analysis: Analysis) -> list[str]:
+    """One line for the analysis and one per channel: samples per kind, processes, regions, categories, variables,
+    variations and estimators."""
+    lines = [f"analysis {analysis.name}: era {analysis.era}, {analysis.lumi_pb:g} pb^-1, signal {analysis.signal}, channels {', '.join(analysis.channels)}"]
+    for name, channel in analysis.channels.items():
+        kinds = ", ".join(f"{kind} {count}" for kind, count in sorted(Counter(s.kind for s in channel.samples).items()))
+        lines.append(
+            f"{name}: {len(channel.samples)} samples ({kinds}), {len(channel.processes)} processes, {len(channel.regions)} regions, "
+            f"{len(channel.categories)} categories, {len(channel.variables)} control variables, {len(channel.variations)} variations, "
+            f"estimators {', '.join(type(e).__name__ for e in channel.estimators) or 'none'}"
+        )
+    return lines
+
+
 def load_analysis(config: RunConfig) -> Analysis:
     """Import `module:function` from the config, build the analysis and validate it."""
     module_name, _, function_name = config.analysis.partition(":")
@@ -85,4 +109,6 @@ def load_analysis(config: RunConfig) -> Analysis:
     build = getattr(importlib.import_module(module_name), function_name)
     analysis = build(config)
     validate(analysis)
+    for line in summary(analysis):
+        logger.info(line)
     return analysis
diff --git a/src/shapesmith/datacards.py b/src/shapesmith/datacards.py
index 571fcac..4be68d8 100644
--- a/src/shapesmith/datacards.py
+++ b/src/shapesmith/datacards.py
@@ -7,6 +7,7 @@
 """
 from __future__ import annotations
 
+import logging
 from dataclasses import dataclass
 from pathlib import Path
 
@@ -14,6 +15,8 @@
 from shapesmith.model import NOMINAL, Analysis, Category, LnN
 from shapesmith.shapes import shape_name, write_shapes
 
+logger = logging.getLogger(__name__)
+
 
 @dataclass
 class Bin:
@@ -56,6 +59,7 @@ def key_of(process: str, variation: str = NOMINAL_VARIATION) -> HistKey:
 
     nominal = {p: hset[key_of(p)] for p in (analysis.signal, *channel.backgrounds()) if key_of(p) in hset}
     if not nominal or key_of(channel.data()) not in hset:
+        logger.warning(f"{channel_name}/{category.name}/{variable}: no {'processes' if not nominal else 'data'}, bin {bin_name(analysis, channel_name, index)} left out")
         return None
     total = None
     for process, h in nominal.items():
@@ -63,9 +67,12 @@ def key_of(process: str, variation: str = NOMINAL_VARIATION) -> HistKey:
             total = h.copy() if total is None else total.add(h)
     edges = rebin_edges(total, min_background) if total is not None else list(hset[key_of(channel.data())].edges)
     result = Bin(channel_name, bin_name(analysis, channel_name, index), hset[key_of(channel.data())].rebin(edges), {}, {})
+    if len(edges) < len(hset[key_of(channel.data())].edges):
+        logger.debug(f"{result.name}: {len(hset[key_of(channel.data())].edges) - 1} -> {len(edges) - 1} bins (background >= {min_background:g} per bin)")
     for process, h in nominal.items():
         rebinned = h.rebin(edges)
         if rebinned.sum() <= 0.0:
+            logger.debug(f"{result.name}: {process} left out, rate {rebinned.sum():.3g}")
             continue
         result.processes[process] = rebinned
         keys = hset.select(channel=channel_name, category=category.name, process=process, region=NOMINAL, variable=variable)
@@ -145,7 +152,10 @@ def write_datacard(hset: HistogramSet, analysis: Analysis, channels: list[str],
             entries += [(b.name, shape_name(process, variation), varied) for variation, varied in b.variations[process].items()]
     write_shapes(output_dir / "common" / f"htt_input_{analysis.era}.root", entries)
     card = output_dir / "combined.txt"
-    card.write_text(card_text(analysis, bins, systematics))
+    text = card_text(analysis, bins, systematics)
+    card.write_text(text)
+    nuisances = sum(1 for line in text.splitlines() if " lnN " in line or " shape " in line)
+    logger.info(f"{output_dir.name}: {len(bins)} bins, {sum(len(b.processes) for b in bins)} process columns, {nuisances} nuisances -> {card}")
     return card
 
 
diff --git a/src/shapesmith/estimates.py b/src/shapesmith/estimates.py
index c192cb5..3977665 100644
--- a/src/shapesmith/estimates.py
+++ b/src/shapesmith/estimates.py
@@ -55,11 +55,14 @@ def estimate_data_minus(hset: HistogramSet, channel: Channel, estimator: DataMin
     added = []
     for category, variable in categories_and_variables(hset, channel, estimator.region):
         present = {k.variation for process in inputs for k in hset.select(channel=channel.name, category=category, process=process, region=estimator.region, variable=variable)}
-        for variation in (NOMINAL_VARIATION, *sorted(present & column_variations)):
+        variations = (NOMINAL_VARIATION, *sorted(present & column_variations))
+        for variation in variations:
             h = data_minus(hset, channel, category, variable, estimator.region, estimator.subtract, variation).scale(estimator.scale)
             key = HistKey(channel.name, category, estimator.output, NOMINAL, variation, variable)
             hset[key] = clip_negative_bins(h) if estimator.clip_negative else h
             added.append(key)
+        nominal = hset[HistKey(channel.name, category, estimator.output, NOMINAL, NOMINAL_VARIATION, variable)]
+        logger.debug(f"{channel.name}/{category}/{variable}: {estimator.output} = data - {' - '.join(estimator.subtract)} in {estimator.region}, yield {nominal.sum():.2f}, {len(variations) - 1} variations")
     return added
 
 
@@ -90,6 +93,7 @@ def estimate_abcd(hset: HistogramSet, channel: Channel, estimator: ABCD) -> list
             factor = 0.0
         else:
             factor = yield_c / yield_d
+        logger.debug(f"{channel.name}/{category}/{variable}: ABCD {estimator.output}, B={b.sum():.2f}, C={yield_c:.2f}, D={yield_d:.2f}, C/D={factor:.4f}")
         key = HistKey(channel.name, category, estimator.output, NOMINAL, NOMINAL_VARIATION, variable)
         hset[key] = clip_negative_bins(b.scale(factor))
         added.append(key)
@@ -142,5 +146,7 @@ def run_estimates(hset: HistogramSet, analysis: Analysis, channels: list[str]) -
     for name in channels:
         channel = analysis.channel(name)
         for estimator in channel.estimators:
-            added += ESTIMATES[type(estimator)](hset, channel, estimator)
+            keys = ESTIMATES[type(estimator)](hset, channel, estimator)
+            logger.info(f"{name}: {type(estimator).__name__} {getattr(estimator, 'output', None) or estimator.name}, {len(keys)} histograms")
+            added += keys
     return added
diff --git a/src/shapesmith/events.py b/src/shapesmith/events.py
index da89cb3..0566b12 100644
--- a/src/shapesmith/events.py
+++ b/src/shapesmith/events.py
@@ -9,6 +9,7 @@
 """
 from __future__ import annotations
 
+import logging
 from dataclasses import dataclass
 from typing import Iterable
 
@@ -20,6 +21,8 @@
 from shapesmith.model import NOMINAL, Analysis, Channel, ColumnVariation, Process, Region, Variation, WeightVariation
 from shapesmith.store import SKIM_COLUMNS, read_skims, schema
 
+logger = logging.getLogger(__name__)
+
 
 @dataclass(frozen=True)
 class Query:
@@ -71,4 +74,6 @@ def load(config: RunConfig, analysis: Analysis, query: Query, columns: Iterable[
     cuts, weights = selection
     frame = read_skims(config.skim_dir, channel.name, nicks, set(SKIM_COLUMNS) | columns_of(cuts + weights) | set(columns))
     selected = frame[mask(frame, cuts)].reset_index(drop=True)
+    variation = f" ({query.variation.name})" if query.variation is not None else ""
+    logger.debug(f"{query.channel}/{query.process}#{query.region}{variation}: {len(selected)} of {len(frame)} events selected")
     return Events(selected, event_weights(selected, weights, lumi(analysis, channel, process)))
diff --git a/src/shapesmith/fill.py b/src/shapesmith/fill.py
index 29b1871..e5cae79 100644
--- a/src/shapesmith/fill.py
+++ b/src/shapesmith/fill.py
@@ -8,12 +8,14 @@
 from __future__ import annotations
 
 import logging
+import time
 from dataclasses import dataclass
 from pathlib import Path
 from typing import NamedTuple
 
 import numpy as np
 
+from shapesmith import logs
 from shapesmith.binning import RECORD, resolve, write_record
 from shapesmith.config import RunConfig
 from shapesmith.events import event_weights, lumi, select
@@ -130,6 +132,7 @@ def fill_booking(config: RunConfig, analysis: Analysis, channel_name: str, booki
     unresolved = [t.variable.name for t in targets_ if isinstance(t.variable.edges, EqualData)]
     if unresolved:
         raise ValueError(f"{channel_name}: the equal-data binning of {', '.join(unresolved)} is not resolved (fill through run_hist)")
+    start = time.monotonic()
     channel = analysis.channel(channel_name)
     process = booking.process
     nicks = [s.nick for s in channel.samples_of(process.group)]
@@ -152,6 +155,10 @@ def fill_booking(config: RunConfig, analysis: Analysis, channel_name: str, booki
             keep = in_category & finite
             key = HistKey(channel_name, target.category, process.name, f.region, f.variation.name if f.variation else NOMINAL_VARIATION, target.variable.name)
             result[key] = Histogram.fill(target.variable.edges, values[keep], w[keep])
+    logger.debug(
+        f"{channel_name}/{process.name}: {len(result)} histograms from {len(frame)} events of {len(nicks)} samples "
+        f"({len(fills)} fills in {len(booking.regions)} regions), {logs.duration(time.monotonic() - start)}"
+    )
     return result
 
 
@@ -170,23 +177,27 @@ def run_hist(config: RunConfig, analysis: Analysis, channels: list[str] | None,
     for channel_name in channels or config.channels:
         channel = analysis.channel(channel_name)
         targets_, record[channel_name] = resolve(config, analysis, channel_name, targets(channel, control, variables))
-        for booking in bookings(channel, systematics, regions):
-            if processes and booking.process.name not in processes:
-                continue
+        chosen = [b for b in bookings(channel, systematics, regions) if not processes or b.process.name in processes]
+        for booking in chosen:
             jobs.append((config, analysis, channel_name, booking, targets_))
             scopes |= {(channel_name, booking.process.name, region, t.category, t.variable.name) for region in booking.regions for t in targets_}
+        logger.info(
+            f"{channel_name}: {len(chosen)} processes, {len(targets_)} {'control variables' if control else 'categories'}, "
+            f"regions {', '.join(dict.fromkeys(r for b in chosen for r in b.regions))}, {'with' if systematics else 'without'} systematics"
+        )
     hset = HistogramSet()
     if output.exists():
         previous = HistogramSet.load(output)
         hset.update((key, h) for key, h in previous.items() if (key.channel, key.process, key.region, key.category, key.variable) not in scopes)
         logger.info(f"{output}: filling {len(scopes)} histogram scopes, keeping {len(hset)} other histograms")
     logger.info(f"filling {len(jobs)} process bookings with {config.workers} workers")
-    for job, result, error in run_jobs(_fill_job, jobs, config.workers):
+    for job, result, error in run_jobs(_fill_job, jobs, config.workers, label="hist"):
         if error is not None:
             raise RuntimeError(f"{job[2]}/{job[3].process.name}: {error}") from error
         hset.update(result)
     hset.save(output)
+    logger.info(f"{len(hset)} histograms written to {output}")
     record = {channel: entries for channel, entries in record.items() if entries}
     if record:
-        write_record(output.parent / RECORD, record)
+        logger.info(f"equal-data bin edges recorded in {write_record(output.parent / RECORD, record)}")
     return hset
diff --git a/src/shapesmith/histogram.py b/src/shapesmith/histogram.py
index 8daff12..c064305 100644
--- a/src/shapesmith/histogram.py
+++ b/src/shapesmith/histogram.py
@@ -8,6 +8,7 @@
 from __future__ import annotations
 
 import json
+import logging
 from dataclasses import asdict, dataclass
 from pathlib import Path
 
@@ -15,6 +16,8 @@
 import uproot
 from uproot.writing.identify import to_TAxis, to_TH1x
 
+logger = logging.getLogger(__name__)
+
 NOMINAL_VARIATION = "Nominal"
 INCLUSIVE = "inclusive"  # the category of control-variable histograms
 TEMPLATE_SEPARATOR = "@"
@@ -161,6 +164,7 @@ def save(self, path: Path) -> None:
                 f[key.path] = h.to_root(key.object_name)
         index = [{**asdict(key), "sum": h.sum(), "bins": len(h.values)} for key, h in items]
         path.with_suffix(".json").write_text(json.dumps(index, indent=1))
+        logger.debug(f"{len(items)} histograms saved to {path} (+ {path.with_suffix('.json').name})")
 
     @classmethod
     def load(cls, path: Path) -> "HistogramSet":
@@ -169,6 +173,7 @@ def load(cls, path: Path) -> "HistogramSet":
             for name, classname in f.classnames(recursive=True, cycle=False).items():
                 if classname.startswith("TH1"):
                     hset[HistKey.parse(name)] = Histogram.from_root(f[name])
+        logger.debug(f"{len(hset)} histograms loaded from {path}")
         return hset
 
 
diff --git a/src/shapesmith/limits.py b/src/shapesmith/limits.py
index 8f9d93b..84f9837 100644
--- a/src/shapesmith/limits.py
+++ b/src/shapesmith/limits.py
@@ -2,6 +2,7 @@
 from __future__ import annotations
 
 import json
+import logging
 from pathlib import Path
 
 import matplotlib
@@ -16,6 +17,8 @@
 QUANTILES = {0.025: "exp_m2", 0.16: "exp_m1", 0.5: "exp_median", 0.84: "exp_p1", 0.975: "exp_p2", -1.0: "observed"}
 COMMON = "-m {mass} --setParameterRanges r=-40,40 -t -1"
 
+logger = logging.getLogger(__name__)
+
 
 def limit_commands(mass: str = "125") -> list[str]:
     common = COMMON.format(mass=mass)
@@ -39,17 +42,21 @@ def collect(card_dir: Path, mass: str = "125") -> dict:
     card_dir = Path(card_dir)
     result: dict = {}
     limit_file = card_dir / f"higgsCombine.Limit.AsymptoticLimits.mH{mass}.root"
+    significance_file = card_dir / f"higgsCombine.Significance.Significance.mH{mass}.root"
+    fit_file = card_dir / f"higgsCombine.Fit.MultiDimFit.mH{mass}.root"
+    for path in (limit_file, significance_file, fit_file):
+        if not path.exists():
+            logger.warning(f"{card_dir.name}: {path.name} missing")
     if limit_file.exists():
         result["limit"] = {QUANTILES[q]: limit for q, limit, _ in _read(limit_file) if q in QUANTILES}
-    significance_file = card_dir / f"higgsCombine.Significance.Significance.mH{mass}.root"
     if significance_file.exists():
         result["significance"] = _read(significance_file)[0][1]
-    fit_file = card_dir / f"higgsCombine.Fit.MultiDimFit.mH{mass}.root"
     if fit_file.exists():
         entries = _read(fit_file)
         best = [r for q, _, r in entries if q == -1.0]
         others = [r for q, _, r in entries if q != -1.0]
         result["r"] = {"best": best[0] if best else None, "low": min(others) if others else None, "high": max(others) if others else None}
+    logger.debug(f"{card_dir.name}: {json.dumps(result)}")
     return result
 
 
@@ -86,6 +93,7 @@ def write_summary(results: dict[str, dict], output_dir: Path) -> tuple[Path, Pat
     fig.tight_layout()
     fig.savefig(pdf_path)
     plt.close(fig)
+    logger.info(f"summary of {len(results)} final states: {json_path}, {md_path.name}, {pdf_path.name}")
     return json_path, md_path, pdf_path
 
 
@@ -95,6 +103,7 @@ def run_limits(config: RunConfig, datacard_dir: Path, final_states: list[str], s
     for name in final_states:
         card_dir = datacard_dir / name
         if not skip_combine:
+            logger.info(f"{name}: combine in {card_dir}")
             cmssw.run(limit_commands(), config.combine, card_dir)
         results[name] = collect(card_dir)
     write_summary(results, datacard_dir)
diff --git a/src/shapesmith/logs.py b/src/shapesmith/logs.py
new file mode 100644
index 0000000..9ce0f58
--- /dev/null
+++ b/src/shapesmith/logs.py
@@ -0,0 +1,131 @@
+"""Logging: the console (rich, on stderr) and a log file per command, one level for ShapeSmith and the analysis.
+
+The level (`log_level` of the run configuration) applies to the logger `shapesmith` and to those of the analysis
+packages; every other library (matplotlib, numexpr, uproot, ...) logs from WARNING on, so that DEBUG stays readable.
+Python warnings are logged as well. The console shows the messages, the log file also the tracebacks, the time, the
+process and the logger. Worker processes of `parallel.run_jobs` hand their records to the main process through a
+queue (`worker_init`, `forwarding`), so their messages reach the same handlers.
+"""
+from __future__ import annotations
+
+import logging
+import logging.handlers
+import multiprocessing
+from contextlib import contextmanager
+from datetime import datetime
+from pathlib import Path
+from typing import Iterable, Iterator
+
+from rich.console import Console
+from rich.logging import RichHandler
+
+LEVELS = ("DEBUG", "INFO", "WARNING", "ERROR")
+FILE_FORMAT = "%(asctime)s %(levelname)-7s %(processName)s %(name)s: %(message)s"
+CORE = "shapesmith"
+OTHERS = logging.WARNING  # the level of every other logger
+
+_installed: list[logging.Handler] = []  # the handlers of configure(), removed by reset()
+_levels: dict[str, int] = {}  # the package levels of configure(), handed to the workers
+
+
+class _MessageOnly(logging.Formatter):
+    """The message without the traceback: on the console, typer shows the traceback of a failing command."""
+
+    def format(self, record: logging.LogRecord) -> str:
+        return record.getMessage()
+
+
+def configure(level: str = "INFO", log_file: Path | None = None, packages: Iterable[str] = ()) -> Path | None:
+    """Log to the console and, if given, to `log_file`, at `level` for `shapesmith` and `packages`; replaces the
+    handlers of an earlier call."""
+    reset()
+    root = logging.getLogger()
+    console = RichHandler(console=Console(stderr=True), show_path=False, markup=False)
+    console.setFormatter(_MessageOnly())
+    _installed.append(console)
+    if log_file is not None:
+        log_file = Path(log_file)
+        log_file.parent.mkdir(parents=True, exist_ok=True)
+        file = logging.FileHandler(log_file, encoding="utf-8")
+        file.setFormatter(logging.Formatter(FILE_FORMAT))
+        _installed.append(file)
+    for handler in _installed:
+        root.addHandler(handler)
+    root.setLevel(OTHERS)
+    for name in dict.fromkeys((CORE, *packages)):
+        logging.getLogger(name).setLevel(level)
+        _levels[name] = logging.getLogger(name).level
+    logging.captureWarnings(True)
+    return log_file
+
+
+def reset() -> None:
+    """Remove and close the handlers of configure() (the tests configure more than once per process)."""
+    root = logging.getLogger()
+    for handler in _installed:
+        root.removeHandler(handler)
+        handler.close()
+    _installed.clear()
+    _levels.clear()
+    logging.captureWarnings(False)
+
+
+def active() -> bool:
+    return bool(_installed)
+
+
+def log_path(directory: Path, command: str) -> Path:
+    """`<directory>/logs/<command>_<YYYYmmdd-HHMMSS>.log`."""
+    return Path(directory) / "logs" / f"{command}_{datetime.now():%Y%m%d-%H%M%S}.log"
+
+
+def levels() -> dict[str, int]:
+    """The package levels for the workers: those of configure(), else the effective level of `shapesmith`."""
+    return dict(_levels) or {CORE: logging.getLogger(CORE).getEffectiveLevel()}
+
+
+def worker_init(queue, package_levels: dict[str, int]) -> None:
+    """Initializer of a worker process: every record goes to the queue, with the package levels of the main process."""
+    root = logging.getLogger()
+    for handler in list(root.handlers):
+        root.removeHandler(handler)
+    root.addHandler(logging.handlers.QueueHandler(queue))
+    root.setLevel(OTHERS)
+    for name, level in package_levels.items():
+        logging.getLogger(name).setLevel(level)
+    logging.captureWarnings(True)
+
+
+class _ToLoggers(logging.Handler):
+    """Hand a record of a worker to the logger of its name in the main process, and so to its handlers."""
+
+    def emit(self, record: logging.LogRecord) -> None:
+        logging.getLogger(record.name).handle(record)
+
+
+@contextmanager
+def forwarding(context: multiprocessing.context.BaseContext) -> Iterator[tuple]:
+    """The initializer and its arguments for a process pool of `context`, while a listener forwards the records."""
+    queue = context.Queue()
+    listener = logging.handlers.QueueListener(queue, _ToLoggers())
+    listener.start()
+    try:
+        yield worker_init, (queue, levels())
+    finally:
+        listener.stop()
+        queue.close()
+
+
+def progress_steps(total: int, parts: int = 10) -> list[int]:
+    """The numbers of finished jobs at which the progress is logged: about every tenth, and the last."""
+    return sorted({max(1, round(total * k / parts)) for k in range(1, parts + 1)}) if total else []
+
+
+def duration(seconds: float) -> str:
+    if seconds < 60:
+        return f"{seconds:.1f} s"
+    minutes, seconds = divmod(round(seconds), 60)
+    if minutes < 60:
+        return f"{minutes} min {seconds:02d} s"
+    hours, minutes = divmod(minutes, 60)
+    return f"{hours} h {minutes:02d} min"
diff --git a/src/shapesmith/measurements/__init__.py b/src/shapesmith/measurements/__init__.py
index 07851b5..ab24ee9 100644
--- a/src/shapesmith/measurements/__init__.py
+++ b/src/shapesmith/measurements/__init__.py
@@ -6,6 +6,7 @@
 """
 from __future__ import annotations
 
+import logging
 from dataclasses import dataclass
 from pathlib import Path
 from typing import Iterable, Protocol
@@ -15,6 +16,8 @@
 from shapesmith.model import Analysis
 from shapesmith.provenance import record, write_versions
 
+logger = logging.getLogger(__name__)
+
 
 @dataclass(frozen=True)
 class MeasureContext:
@@ -45,6 +48,8 @@ def run_measure(config: RunConfig, analysis: Analysis, channels: list[str], sugg
     if analysis.measurement is None:
         raise ValueError(f"analysis {analysis.name} defines no measurement")
     output = Path(config.output_dir) / analysis.measurement.name / analysis.era
+    mode = "suggest-binning" if suggest_binning else "merge" if merge else "measure"
+    logger.info(f"measurement {analysis.measurement.name} ({mode}) of channels {', '.join(channels)} -> {output}")
     if not suggest_binning:
         write_versions(config, analysis.name, output)
     analysis.measurement.run(MeasureContext(config, analysis, channels, output, suggest_binning, merge))
diff --git a/src/shapesmith/measurements/fake_factors/binning.py b/src/shapesmith/measurements/fake_factors/binning.py
index e3be7c2..574c646 100644
--- a/src/shapesmith/measurements/fake_factors/binning.py
+++ b/src/shapesmith/measurements/fake_factors/binning.py
@@ -8,6 +8,7 @@
 """
 from __future__ import annotations
 
+import logging
 from typing import TYPE_CHECKING, Iterator
 
 import numpy as np
@@ -16,6 +17,8 @@
     from shapesmith.measurements.fake_factors.measure import EventSource
     from shapesmith.measurements.fake_factors.model import Binned, Equipopulated, FakeFactorMeasurement, Split
 
+logger = logging.getLogger(__name__)
+
 
 def equipopulated(x: np.ndarray, weights: np.ndarray, options: Equipopulated, n_bins: int) -> list[float]:
     inside = (x >= options.low) & (x <= options.high)
@@ -55,6 +58,6 @@ def suggestions(measurement: FakeFactorMeasurement, source: EventSource) -> list
     return result
 
 
-def print_suggestions(measurement: FakeFactorMeasurement, source: EventSource) -> None:
+def log_suggestions(measurement: FakeFactorMeasurement, source: EventSource) -> None:
     for name, category, suggested, current in suggestions(measurement, source):
-        print(f"{source.channel} {name} category {category}: suggested {suggested}, current {current}")
+        logger.info(f"{source.channel} {name} category {category}: suggested {suggested}, current {current}")
diff --git a/src/shapesmith/measurements/fake_factors/measure.py b/src/shapesmith/measurements/fake_factors/measure.py
index 9cea16e..8513958 100644
--- a/src/shapesmith/measurements/fake_factors/measure.py
+++ b/src/shapesmith/measurements/fake_factors/measure.py
@@ -24,7 +24,7 @@
 from shapesmith import payloads
 from shapesmith.events import Events, Query
 from shapesmith.measurements.fake_factors import payload, plots
-from shapesmith.measurements.fake_factors.binning import print_suggestions
+from shapesmith.measurements.fake_factors.binning import log_suggestions
 from shapesmith.measurements.fake_factors.fit import Fitted, fit
 from shapesmith.measurements.fake_factors.hist import Hist
 
@@ -59,6 +59,7 @@ def __init__(self, context: MeasureContext, channel: str, columns: set[str]):
     def events(self, process: str, region: str) -> Events:
         if (process, region) not in self._events:
             self._events[process, region] = self._context.events(Query(self.channel, process, region), self._columns)
+            logger.debug(f"{self.channel}: {process} in {region}, {len(self._events[process, region].weights)} events, sum of weights {self._events[process, region].weights.sum():.4g}")
         return self._events[process, region]
 
     def hist(self, process: str, region: str, variable: str, edges, split: Split, category: int, weights: np.ndarray | None = None) -> Hist:
@@ -127,6 +128,7 @@ def process_fake_factors(self, name: str, process: ProcessFF, sr_like: str, ar_l
         MC fake factor."""
         scale = None if process.scale is None else self.data_scale(process)
         if scale is not None:
+            logger.info(f"{self.source.channel} {stage}{name}: data/MC scale factor {scale:.4f}")
             self.record[f"{stage}{name}_data_scale"] = [{"factor": scale}]
 
         def measure(category: int) -> Ratio:
@@ -238,6 +240,12 @@ def fitted(self, name: str, binned: Binned, split: Split, measure: Callable[[int
         for category in split.categories:
             ratio = measure(category)
             fitted = fit(ratio.nominal, binned.fits[category], ratio.mc_shifted, self.measurement.stat_sigma)
+            logger.debug(
+                f"{self.source.channel} {name} category {category} ({split.variable} {split.edges[category]}-{split.edges[category + 1]}): {binned.variable} fit {binned.fits[category]}, "
+                f"ratio {', '.join(f'{v:.3g}' for v in ratio.nominal.values)}" + (f", p-value {fitted.p_value:.3g}" if fitted.p_value is not None else "")
+            )
+            if fitted.reset:
+                logger.info(f"{self.source.channel} {name} category {category}: compatible with 1 (p-value {fitted.p_value:.3g}), set to 1")
             self.record.setdefault(name, []).append(_record(binned.variable, split, category, ratio, fitted))
             result.append(fitted)
         return result
@@ -247,11 +255,12 @@ def run(measurement: FakeFactorMeasurement, context: MeasureContext) -> None:
     for channel in context.channels:
         source = EventSource(context, channel, measurement.columns(channel))
         if context.suggest_binning:
-            print_suggestions(measurement, source)
+            log_suggestions(measurement, source)
             continue
-        logger.info(f"{channel}: measuring the fake factors")
+        logger.info(f"{channel}: measuring the fake factors of {len(measurement.legs[channel])} legs")
         result = ChannelMeasurement(measurement, source)
         result.measure(measurement.legs[channel])
+        logger.info(f"{channel}: {len(result.fake_factors)} fake factors and fractions, {len(result.corrections)} corrections, {len(result.intermediate)} intermediate corrections")
         provenance, output = context.provenance(channel=channel), context.output
         payloads.write(payloads.correction_set(result.fake_factors, provenance), output / f"fake_factors_{channel}.json.gz")
         payloads.write(payloads.correction_set(result.corrections, provenance, result.compounds), output / f"FF_corrections_{channel}.json.gz")
diff --git a/src/shapesmith/measurements/tau_id_es/__init__.py b/src/shapesmith/measurements/tau_id_es/__init__.py
index 4401c80..00eb517 100644
--- a/src/shapesmith/measurements/tau_id_es/__init__.py
+++ b/src/shapesmith/measurements/tau_id_es/__init__.py
@@ -98,6 +98,7 @@ def fit_or_failure(combine_config: CombineConfig | None, category: str, synced:
     try:
         return fit_category(combine_config, category, synced, directory, era, grid)
     except (subprocess.CalledProcessError, subprocess.TimeoutExpired) as error:
+        logger.error(f"{category}: CMSSW step failed: {error}")
         return {"sf": None, "es": None, "problems": [f"CMSSW step failed: {error}"]}
 
 
@@ -114,17 +115,21 @@ def measure(measurement: TauIdEsMeasurement, context: MeasureContext) -> Path:
     check(analysis, measurement)
     directory = context.output / f"{measurement.vsjet_wp}_{measurement.vsele_wp}"
     shapes = directory / "shapes.root"
+    logger.info(f"tau_id_es {measurement.vsjet_wp}/{measurement.vsele_wp}: histograms of {CHANNEL} and {CONTROL_CHANNEL} with {len(measurement.grid)} ES shifts -> {directory}")
     hset = run_hist(config, analysis, [CHANNEL, CONTROL_CHANNEL], False, None, True, None, shapes)
     run_estimates(hset, analysis, [CHANNEL, CONTROL_CHANNEL])
     hset.save(shapes)
     synced = write_synced(hset, analysis, directory / "synced")
+    logger.info(f"synced shapes: {synced}")
     categories = [c.name for c in analysis.channel(CHANNEL).categories]
     for category in categories:
         plots.plot_control(hset, analysis, category, directory / "plots" / "control" / category)
+    logger.info(f"fitting {len(categories)} categories, {max(1, config.workers)} at a time: {', '.join(categories)}")
     with ThreadPoolExecutor(max_workers=max(1, config.workers)) as pool:  # each category runs in its own CMSSW process
         fits = list(pool.map(lambda category: fit_or_failure(config.combine, category, synced, directory, analysis.era, measurement.grid), categories))
     for category, fit in zip(categories, fits):
         if fit["sf"] is not None:
+            logger.info(f"{category}: SF {_interval(fit['sf'])}, ES {_interval(fit['es'])}" + (f", problems: {'; '.join(fit['problems'])}" if fit["problems"] else ""))
             plot_category(directory, category, fit)
     record = {
         "vsjet_wp": measurement.vsjet_wp,
@@ -139,9 +144,14 @@ def measure(measurement: TauIdEsMeasurement, context: MeasureContext) -> Path:
     return directory
 
 
+def _interval(interval: dict) -> str:
+    return " ".join(f"{name}={value:.4g}" for name, value in interval.items() if isinstance(value, (int, float)))
+
+
 def merge(context: MeasureContext) -> Path:
     """The payload of every working-point combination measured below the output directory."""
     records = [json.loads(path.read_text()) for path in sorted(Path(context.output).glob("*/results.json"))]
+    logger.info(f"tau_id_es: merging {len(records)} working-point combinations from {context.output}")
     if not records:
         raise ValueError(f"tau_id_es: no results.json below {context.output}; run the measurement per working point first")
     problems = [f"{r['vsjet_wp']}_{r['vsele_wp']}/{category}: {problem}" for r in records for category, fit in r["categories"].items() for problem in fit["problems"]]
diff --git a/src/shapesmith/ml_export.py b/src/shapesmith/ml_export.py
index bbf43ea..d3d0da0 100644
--- a/src/shapesmith/ml_export.py
+++ b/src/shapesmith/ml_export.py
@@ -100,5 +100,7 @@ def run_ml_export(config: RunConfig, analysis: Analysis, channels: list[str]) ->
             combined[tuple_column("Nominal", "class_weight")] = class_weights(combined[tuple_column("Nominal", "weight")].to_numpy().astype(np.float64), labels).astype(np.float32)
             path = directory / f"{fold}.feather"
             feather.write_feather(pa.Table.from_pandas(combined), path)  # pyarrow keeps the MultiIndex columns in its pandas metadata
+            logger.debug(f"{channel}/{fold}: {len(combined)} events -> {path}")
             written.append(path)
+        logger.info(f"{channel}: {len(parts)} fold files in {directory}")
     return written
diff --git a/src/shapesmith/ntuples.py b/src/shapesmith/ntuples.py
index 7426e5d..bd5e422 100644
--- a/src/shapesmith/ntuples.py
+++ b/src/shapesmith/ntuples.py
@@ -2,6 +2,7 @@
 from __future__ import annotations
 
 import json
+import logging
 import os
 from dataclasses import dataclass
 
@@ -12,6 +13,8 @@
 
 TREE = "ntuple"
 
+logger = logging.getLogger(__name__)
+
 
 @dataclass(frozen=True)
 class NtupleFile:
@@ -70,6 +73,7 @@ def discover(ntuples: NtupleConfig, era: str, nick: str, channel: str, kind: str
     missing = [f"{directory}/{name}" for directory, names in zip(friend_dirs, friend_listings) for name in basenames if name not in names]
     if missing:
         raise FileNotFoundError("missing friend files:\n" + "\n".join(missing))
+    logger.debug(f"{channel}/{nick}: {len(basenames)} files in {join_url(ntuples.server, main_dir)}" + (f", {len(friend_dirs)} friends each" if friend_dirs else ""))
     return [
         NtupleFile(
             path=join_url(ntuples.server, f"{main_dir}/{name}"),
diff --git a/src/shapesmith/parallel.py b/src/shapesmith/parallel.py
index 0b7ae0b..bca69dd 100644
--- a/src/shapesmith/parallel.py
+++ b/src/shapesmith/parallel.py
@@ -2,31 +2,58 @@
 
 `spawn` instead of `fork`: uproot, pyarrow and XRootD start threads, and forking a multi-threaded
 process can deadlock the children. With one worker everything runs in the calling process, which
-keeps tracebacks and debuggers usable (`--workers 1` / `workers: 1`).
+keeps tracebacks and debuggers usable (`--workers 1` / `workers: 1`). The workers log through the
+main process (logs.forwarding); the progress is logged about every tenth of the jobs.
 """
 from __future__ import annotations
 
+import logging
 import multiprocessing
+import time
 from concurrent.futures import ProcessPoolExecutor, as_completed
 from typing import Callable, Iterable, Iterator, TypeVar
 
+from shapesmith import logs
+
 T = TypeVar("T")
 
+logger = logging.getLogger(__name__)
+
 
-def run_jobs(function: Callable[[tuple], T], jobs: Iterable[tuple], workers: int) -> Iterator[tuple[tuple, T | None, Exception | None]]:
+def run_jobs(function: Callable[[tuple], T], jobs: Iterable[tuple], workers: int, label: str = "jobs") -> Iterator[tuple[tuple, T | None, Exception | None]]:
     """Yield (job, result, error) for every job; errors are handed back instead of raised so that callers can report all of them."""
     jobs = list(jobs)
+    steps, start, done = set(logs.progress_steps(len(jobs))), time.monotonic(), 0
+
+    def progress() -> None:
+        nonlocal done
+        done += 1
+        if done in steps:
+            logger.info(f"{label}: {done}/{len(jobs)} jobs done ({100 * done / len(jobs):.0f} %), {logs.duration(time.monotonic() - start)}")
+
     if workers <= 1 or len(jobs) <= 1:
+        logger.debug(f"{label}: {len(jobs)} jobs inline")
         for job in jobs:
             try:
-                yield job, function(job), None
+                result = function(job)
             except Exception as error:  # reported by the caller
+                progress()
                 yield job, None, error
+                continue
+            progress()
+            yield job, result, None
         return
-    with ProcessPoolExecutor(max_workers=min(workers, len(jobs)), mp_context=multiprocessing.get_context("spawn")) as pool:
-        futures = {pool.submit(function, job): job for job in jobs}
-        for future in as_completed(futures):
-            try:
-                yield futures[future], future.result(), None
-            except Exception as error:
-                yield futures[future], None, error
+    context = multiprocessing.get_context("spawn")
+    logger.debug(f"{label}: {len(jobs)} jobs in {min(workers, len(jobs))} worker processes")
+    with logs.forwarding(context) as (initializer, initargs):
+        with ProcessPoolExecutor(max_workers=min(workers, len(jobs)), mp_context=context, initializer=initializer, initargs=initargs) as pool:
+            futures = {pool.submit(function, job): job for job in jobs}
+            for future in as_completed(futures):
+                try:
+                    result = future.result()
+                except Exception as error:
+                    progress()
+                    yield futures[future], None, error
+                    continue
+                progress()
+                yield futures[future], result, None
diff --git a/src/shapesmith/plotting/stack.py b/src/shapesmith/plotting/stack.py
index d55d39f..00be1e7 100644
--- a/src/shapesmith/plotting/stack.py
+++ b/src/shapesmith/plotting/stack.py
@@ -118,6 +118,7 @@ def plot_stack(hset: HistogramSet, analysis: Analysis, channel: str, category: s
         output_dir = output_dir / region
     output_dir.mkdir(parents=True, exist_ok=True)
     paths = [output_dir / f"{category}_{variable}.{ext}" for ext in ("pdf", "png")]
+    logger.debug(f"{channel}/{category}/{variable}#{region}: {len(groups)} background groups, signal x {scale:g} -> {output_dir}/{category}_{variable}.pdf/png")
     for path in paths:
         fig.savefig(path, bbox_inches="tight")
     plt.close(fig)
@@ -128,6 +129,7 @@ def run_plot(hset: HistogramSet, analysis: Analysis, channels: list[str], contro
     written = []
     for channel in channels:
         analysis.channel(channel).region(region)
+        before = len(written)
         if control:
             for variable in variables or list(analysis.channel(channel).variables):
                 written += plot_stack(hset, analysis, channel, INCLUSIVE, variable, output_dir, region=region, **options)
@@ -135,4 +137,5 @@ def run_plot(hset: HistogramSet, analysis: Analysis, channels: list[str], contro
             for cat in analysis.channel(channel).categories:
                 if category is None or cat.name == category:
                     written += plot_stack(hset, analysis, channel, cat.name, cat.variable.name, output_dir, region=region, **options)
+        logger.info(f"{channel}: {(len(written) - before) // 2} plots of region {region} in {Path(output_dir) / channel}")
     return written
diff --git a/src/shapesmith/provenance.py b/src/shapesmith/provenance.py
index b934d75..5d58b65 100644
--- a/src/shapesmith/provenance.py
+++ b/src/shapesmith/provenance.py
@@ -4,6 +4,7 @@
 import hashlib
 import importlib
 import json
+import logging
 import subprocess
 from datetime import datetime, timezone
 from pathlib import Path
@@ -11,6 +12,8 @@
 from shapesmith import __version__
 from shapesmith.config import RunConfig
 
+logger = logging.getLogger(__name__)
+
 
 def git_hash(path: Path) -> str | None:
     """HEAD commit of the git checkout containing `path`, or None if there is none (or git is missing)."""
@@ -42,4 +45,5 @@ def write_versions(config: RunConfig, analysis_name: str, directory: Path) -> Pa
     path = Path(directory) / "versions" / f"{analysis_name}.json"
     path.parent.mkdir(parents=True, exist_ok=True)
     path.write_text(json.dumps(record(config), indent=2))
+    logger.debug(f"versions recorded in {path}")
     return path
diff --git a/src/shapesmith/shapes.py b/src/shapesmith/shapes.py
index 2f2b9c7..fe5d897 100644
--- a/src/shapesmith/shapes.py
+++ b/src/shapesmith/shapes.py
@@ -1,6 +1,7 @@
 """combine-style shape files: one writer, and the synced layout <channel>_<category>/{data_obs, <process>, <process>_<variation>}."""
 from __future__ import annotations
 
+import logging
 from pathlib import Path
 from typing import Callable, Iterable
 
@@ -9,14 +10,19 @@
 from shapesmith.histogram import NOMINAL_VARIATION, HistKey, Histogram, HistogramSet, is_part
 from shapesmith.model import NOMINAL, Analysis
 
+logger = logging.getLogger(__name__)
+
 
 def write_shapes(path: Path, entries: Iterable[tuple[str, str, Histogram]]) -> Path:
     """Write (folder, name, histogram) entries as TH1D `folder/name` into a new ROOT file."""
     path = Path(path)
     path.parent.mkdir(parents=True, exist_ok=True)
+    count = 0
     with uproot.recreate(path) as f:
         for folder, name, h in entries:
             f[f"{folder}/{name}"] = h.to_root(name)
+            count += 1
+    logger.debug(f"{count} histograms written to {path}")
     return path
 
 
diff --git a/src/shapesmith/skim.py b/src/shapesmith/skim.py
index 74c4d8d..7033310 100644
--- a/src/shapesmith/skim.py
+++ b/src/shapesmith/skim.py
@@ -14,6 +14,7 @@
 from __future__ import annotations
 
 import logging
+import time
 from concurrent.futures import ThreadPoolExecutor
 from dataclasses import dataclass, field
 from pathlib import Path
@@ -21,12 +22,12 @@
 import numpy as np
 import pandas as pd
 
-from shapesmith import __version__
+from shapesmith import __version__, logs
 from shapesmith.config import RunConfig
 from shapesmith.events import select
 from shapesmith.expressions import columns_in, columns_of, mask, shift
 from shapesmith.model import NOMINAL_REGION, Analysis, Channel, ColumnVariation, Sample, WeightVariation, applies
-from shapesmith.ntuples import NtupleFile, discover, friend_bases, read_ntuple
+from shapesmith.ntuples import NtupleFile, discover, friend_bases, join_url, read_ntuple
 from shapesmith.parallel import run_jobs
 from shapesmith.store import SKIM_COLUMNS, manifest_path, read_manifest, schema, skim_path, write_manifest, write_skim
 
@@ -87,9 +88,11 @@ def check_shifts(channel: Channel, sample: Sample, columns: set[str], branches:
 def skim_one(ntuple: NtupleFile, sample: Sample, channel: Channel, columns: set[str], out_path: Path, optional: set[str] = frozenset()) -> SkimResult:
     """Read one file, keep the events that pass the skim cuts under any variation and the sample cut, add the
     bookkeeping columns and write Parquet. Every column in `columns` must exist; `optional` ones may be absent."""
+    start = time.monotonic()
     variations = column_variations(channel, sample)
     shifted = {column + v.suffix for v in variations if v.suffix for column in columns}
     frame, metadata, branches = read_ntuple(ntuple, columns | shifted, set(optional) | shifted)
+    read = time.monotonic() - start
     check_shifts(channel, sample, columns, branches)
     n_in = len(frame)
     nominal = list(channel.skim.values())
@@ -108,6 +111,10 @@ def skim_one(ntuple: NtupleFile, sample: Sample, channel: Channel, columns: set[
     for column, kind in (("is_data", "data"), ("is_mc", "mc"), ("is_embedding", "embedding")):
         frame[column] = np.full(len(frame), sample.kind == kind, dtype=bool)
     write_skim(frame, out_path)
+    logger.debug(
+        f"{channel.name}/{sample.nick}/{ntuple.basename}: {n_in} -> {len(frame)} events, {len(frame.columns)} columns "
+        f"({len(shifted & set(frame.columns))} shifted), read in {logs.duration(read)}, done in {logs.duration(time.monotonic() - start)}"
+    )
     return SkimResult(sample.nick, ntuple.basename, n_in, len(frame), out_path, metadata)
 
 
@@ -181,8 +188,10 @@ def _plan(config: RunConfig, channel: Channel, sample: Sample, files: list[Ntupl
             return manifest, [], ", ".join(problems)
         completed = {name: record for name, record in previous.get("completed", {}).items() if name in {f.basename for f in files}}
         if all(ntuple.basename in completed and out.exists() for ntuple, out in outputs):
+            logger.debug(f"{channel.name}/{sample.nick}: stored skim of {len(files)} files reused")
             return {**previous, "completed": completed}, [], None
         if all(previous["contract"].get(name) == manifest["contract"].get(name) for name in ROW_FIELDS):
+            logger.debug(f"{channel.name}/{sample.nick}: resuming the stored skim, {len(completed)} of {len(files)} files done")
             manifest.update(metadata=previous.get("metadata", {}), completed=completed)
             columns |= stored_columns
             manifest["contract"]["required_columns"] = sorted(columns)
@@ -190,6 +199,7 @@ def _plan(config: RunConfig, channel: Channel, sample: Sample, files: list[Ntupl
             logger.info(f"{channel.name}/{sample.nick}: the stored skim holds more variations but is incomplete, skimming it again")
             stored_columns = set()
     jobs = [(ntuple, sample, channel, columns, out, stored_columns) for ntuple, out in outputs if not (ntuple.basename in manifest["completed"] and out.exists())]
+    logger.debug(f"{channel.name}/{sample.nick}: {len(jobs)} of {len(files)} files to skim, {len(columns)} columns" + (" (--force)" if force else ""))
     for ntuple, *_ in jobs:
         manifest["completed"].pop(ntuple.basename, None)  # its Parquet is gone: the record no longer counts
     return manifest, jobs, None
@@ -219,7 +229,11 @@ def run_skim(config: RunConfig, analysis: Analysis, channels: list[str] | None =
         chosen = select_samples(channel, samples)
         with ThreadPoolExecutor(max_workers=min(8, max(1, len(chosen)))) as listing:  # directory listings are latency bound
             discovered = list(listing.map(lambda s: discover(config.ntuples, config.era, s.nick, channel_name, s.kind), chosen))
+        logger.info(f"{channel_name}: {len(chosen)} samples, {sum(len(files) for files in discovered)} files in {join_url(config.ntuples.server, config.ntuples.base)}")
+        counts = dict.fromkeys(("reused", "resumed", "new", "empty", "pending"), 0)
         for sample, files in zip(chosen, discovered):
+            if not files:
+                logger.warning(f"{channel_name}/{sample.nick}: no ntuple files")
             manifest, pending, problem = _plan(config, channel, sample, files, force)
             if problem is not None:
                 incompatible.append(f"{channel_name}/{sample.nick}: {problem}")
@@ -227,17 +241,21 @@ def run_skim(config: RunConfig, analysis: Analysis, channels: list[str] | None =
             manifests[channel_name, sample.nick] = manifest
             results += [SkimResult(sample.nick, name, record["n_in"], record["n_out"], skim_path(config.skim_dir, channel_name, sample.nick, name), skipped=True) for name, record in manifest["completed"].items()]
             jobs += pending
+            counts["empty" if not files else "reused" if not pending else "resumed" if manifest["completed"] else "new"] += 1
+            counts["pending"] += len(pending)
+        logger.info(f"{channel_name}: samples {counts['reused']} reused, {counts['resumed']} resumed, {counts['new']} new, {counts['empty']} without files; {counts['pending']} files to skim")
     if samples and not manifests and not incompatible:
         raise ValueError(f"no sample matches {samples}")
     if incompatible:
         raise ValueError("incompatible skims (re-run `shapesmith skim --force --samples ...` for them):\n" + "\n".join(incompatible))
     for channel_name, nick in {(job[2].name, job[1].nick) for job in jobs}:
         _write_manifest(config, manifests[channel_name, nick])  # the contract is recorded before the first file is written
-    logger.info(f"skimming {len(jobs)} files with {config.workers} workers")
+    logger.info(f"skimming {len(jobs)} files with {config.workers} workers into {config.skim_dir}")
     failures = []
-    for job, result, error in run_jobs(_skim_job, jobs, config.workers):
+    for job, result, error in run_jobs(_skim_job, jobs, config.workers, label="skim"):
         if error is not None:  # collect, report all at the end
             failures.append(f"{job[0].path}: {error}")
+            logger.error(f"{job[2].name}/{job[1].nick}/{job[0].basename}: {error}")
             continue
         results.append(result)
         manifest = manifests[job[2].name, result.nick]
@@ -245,6 +263,9 @@ def run_skim(config: RunConfig, analysis: Analysis, channels: list[str] | None =
         if not manifest["metadata"]:
             manifest["metadata"] = result.metadata
         _write_manifest(config, manifest)
+    for channel_name in channels or config.channels:
+        done = [r for r in results if r.path.parent.parent.name == channel_name]  # <skim_dir>/<channel>/<nick>/<file>
+        logger.info(f"{channel_name}: {len(done)} files, {sum(r.n_in for r in done)} -> {sum(r.n_out for r in done)} events")
     if failures:
-        raise RuntimeError("skim failures:\n" + "\n".join(failures))
+        raise RuntimeError(f"{len(failures)} skim failures:\n" + "\n".join(failures))
     return results
diff --git a/src/shapesmith/store.py b/src/shapesmith/store.py
index 4b2e9ee..5a4ec4d 100644
--- a/src/shapesmith/store.py
+++ b/src/shapesmith/store.py
@@ -2,6 +2,7 @@
 from __future__ import annotations
 
 import json
+import logging
 import os
 from pathlib import Path
 from typing import Callable, Iterable
@@ -14,6 +15,8 @@
 
 SKIM_COLUMNS = ("sample_nick", "norm_weight", "is_data", "is_mc", "is_embedding")  # added to every skim
 
+logger = logging.getLogger(__name__)
+
 
 class SkimMissingError(FileNotFoundError):
     pass
@@ -65,6 +68,7 @@ def schema(skim_dir: Path, channel: str, nicks: Iterable[str]) -> set[str]:
 def read_skims(skim_dir: Path, channel: str, nicks: Iterable[str], columns: Iterable[str] | None) -> pd.DataFrame:
     """All skim rows of the given nicks in one DataFrame, in nick and file order (only `columns` if given)."""
     frames = []
+    nicks = list(nicks)
     for nick in nicks:
         dataset = ds.dataset([str(f) for f in _files(skim_dir, channel, nick)], format="parquet")
         if columns is not None:
@@ -73,7 +77,9 @@ def read_skims(skim_dir: Path, channel: str, nicks: Iterable[str], columns: Iter
                 raise SkimMissingError(f"skims of {nick} in {Path(skim_dir) / channel / nick} lack columns {missing} (re-run `shapesmith skim --force` after changing the analysis)")
         table = dataset.to_table(columns=list(columns) if columns is not None else None)
         frames.append(table.to_pandas())
-    return pd.concat(frames, ignore_index=True) if frames else pd.DataFrame()
+    frame = pd.concat(frames, ignore_index=True) if frames else pd.DataFrame()
+    logger.debug(f"{channel}: {len(frame)} skim rows, {len(frame.columns)} columns of {len(nicks)} samples read")
+    return frame
 
 
 def write_manifest(path: Path, manifest: dict) -> None:
diff --git a/tests/test_cli.py b/tests/test_cli.py
index 1087337..fc5d03e 100644
--- a/tests/test_cli.py
+++ b/tests/test_cli.py
@@ -3,6 +3,7 @@
 import uproot
 from typer.testing import CliRunner
 
+from shapesmith import __version__
 from shapesmith.cli import app
 from shapesmith.testing import make_mini_dataset
 
@@ -76,3 +77,26 @@ def test_validate_reports_broken_analysis(tmp_path):
     config.write_text(config.read_text().replace("tests.mini_analysis:build", "tests.mini_analysis:no_such"))
     result = runner.invoke(app, ["validate", "-c", str(config)])
     assert result.exit_code != 0
+
+
+def test_commands_write_a_log_file(tmp_path):
+    config = _write_config(tmp_path)
+    _run("validate", "-c", str(config))
+    assert not (tmp_path / "out" / "logs").exists()  # validate only reports
+    _run("skim", "-c", str(config), "-s", "log_level=debug")
+    [log] = (tmp_path / "out" / "logs").glob("skim_*.log")
+    text = log.read_text()
+    assert f"shapesmith.cli: ShapeSmith {__version__} skim: {config}, overrides log_level=debug" in text
+    assert "log level DEBUG" in text and "shapesmith.skim: mt: samples 0 reused, 0 resumed" in text
+    assert "shapesmith.parallel: skim: " in text and "jobs done (100 %)" in text
+    assert any(" DEBUG " in line and "SpawnProcess" in line and "shapesmith.skim: mt/" in line for line in text.splitlines())  # from the workers
+    assert "shapesmith.cli: done in" in text
+
+
+def test_a_failing_command_logs_its_traceback(tmp_path):
+    config = _write_config(tmp_path)
+    result = runner.invoke(app, ["hist", "-c", str(config)])  # no skims yet
+    assert result.exit_code != 0
+    [log] = (tmp_path / "out" / "logs").glob("hist_*.log")
+    text = log.read_text()
+    assert "ERROR   MainProcess shapesmith.cli: failed after" in text and "SkimMissingError" in text and "Traceback" in text
diff --git a/tests/test_fake_factors.py b/tests/test_fake_factors.py
index 444754c..53c8d90 100644
--- a/tests/test_fake_factors.py
+++ b/tests/test_fake_factors.py
@@ -2,6 +2,7 @@
 the chaining of the corrections and the payload conventions of the CROWN fake-factor friend."""
 import gzip
 import json
+import logging
 
 import correctionlib
 import numpy as np
@@ -155,12 +156,12 @@ def test_the_payloads_follow_the_friend_conventions(measured):
     assert (output / "FF_for_DRtoSR_mt.json.gz").exists() and list((output / "plots" / "mt").glob("*.png"))
 
 
-def test_suggest_binning_writes_nothing(measured, capsys, tmp_path):
+def test_suggest_binning_writes_nothing(measured, caplog, tmp_path):
     config, analysis, _, _, _ = measured
     config = config.model_copy(update={"output_dir": tmp_path / "suggest"})
-    run_measure(config, analysis, ["mt"], suggest_binning=True)
-    printed = capsys.readouterr().out
-    assert "mt QCD_fake_factors category 0: suggested [30.0," in printed and not (tmp_path / "suggest").exists()
+    with caplog.at_level(logging.INFO, logger="shapesmith"):
+        run_measure(config, analysis, ["mt"], suggest_binning=True)
+    assert "mt QCD_fake_factors category 0: suggested [30.0," in caplog.text and not (tmp_path / "suggest").exists()
 
 
 def test_merge_is_refused(measured, tmp_path):
diff --git a/tests/test_limits.py b/tests/test_limits.py
index fd69f8f..85f09c4 100644
--- a/tests/test_limits.py
+++ b/tests/test_limits.py
@@ -1,3 +1,4 @@
+import logging
 import subprocess
 import time
 from pathlib import Path
@@ -45,6 +46,17 @@ def test_cmssw_run_stops_the_commands_at_the_timeout(tmp_path, monkeypatch):
         cmssw.run(["exit 3"], CombineConfig(cmssw_dir="/unused"), tmp_path)
 
 
+def test_cmssw_run_logs_the_output(tmp_path, monkeypatch, caplog):
+    monkeypatch.setattr(cmssw, "script", lambda combine, cwd, commands: "; ".join([f"cd {cwd}", *commands]))
+    with caplog.at_level(logging.DEBUG, logger="shapesmith"):
+        cmssw.run(["echo hello"], CombineConfig(cmssw_dir="/unused"), tmp_path)
+        with pytest.raises(subprocess.CalledProcessError):
+            cmssw.run(["echo broken >&2", "exit 3"], CombineConfig(cmssw_dir="/unused"), tmp_path)
+    assert "cmssw: hello" in caplog.text
+    [error] = [record for record in caplog.records if record.levelno == logging.ERROR]
+    assert "exit code 3" in error.getMessage() and error.getMessage().endswith("last output:\nbroken")
+
+
 def test_cmssw_run_needs_a_combine_configuration(tmp_path):
     with pytest.raises(ValueError, match="needs `combine`"):
         cmssw.run(["true"], None, tmp_path)
diff --git a/tests/test_logs.py b/tests/test_logs.py
new file mode 100644
index 0000000..04865f1
--- /dev/null
+++ b/tests/test_logs.py
@@ -0,0 +1,99 @@
+import logging
+
+import pytest
+from pydantic import ValidationError
+
+from shapesmith import logs
+from shapesmith.config import NtupleConfig, RunConfig
+from shapesmith.parallel import run_jobs
+
+
+@pytest.fixture
+def restore_logging():
+    """configure() changes the root logger and the package levels; restore them afterwards."""
+    root = logging.getLogger()
+    handlers, level = list(root.handlers), root.level
+    levels = {name: logging.getLogger(name).level for name in ("shapesmith", "analysis_pkg", "matplotlib")}
+    yield
+    logs.reset()
+    root.handlers[:] = handlers
+    root.setLevel(level)
+    for name, value in levels.items():
+        logging.getLogger(name).setLevel(value)
+
+
+def _config(tmp_path, **changes):
+    return RunConfig(analysis="tests.mini_analysis:build", era="2018", channels=["mt"], ntuples=NtupleConfig(base=str(tmp_path)), skim_dir=tmp_path / "skims", output_dir=tmp_path / "out", **changes)
+
+
+def test_log_level_of_the_config(tmp_path):
+    assert _config(tmp_path).log_level == "INFO"
+    assert _config(tmp_path, log_level="debug").log_level == "DEBUG"
+    with pytest.raises(ValidationError):
+        _config(tmp_path, log_level="verbose")
+
+
+def test_configure_sets_the_package_levels_and_writes_the_file(tmp_path, restore_logging):
+    path = logs.configure("DEBUG", tmp_path / "logs" / "run.log", ["analysis_pkg"])
+    logging.getLogger("shapesmith.skim").debug("from the core")
+    logging.getLogger("analysis_pkg.samples").debug("from the analysis")
+    logging.getLogger("matplotlib.font_manager").info("third-party chatter")
+    logging.getLogger("matplotlib").warning("third-party warning")
+    logs.reset()
+    text = path.read_text()
+    assert "DEBUG   MainProcess shapesmith.skim: from the core" in text
+    assert "analysis_pkg.samples: from the analysis" in text
+    assert "third-party chatter" not in text and "third-party warning" in text
+
+
+def test_level_filters_the_shapesmith_messages(tmp_path, restore_logging):
+    path = logs.configure("WARNING", tmp_path / "run.log")
+    logging.getLogger("shapesmith.fill").info("hidden")
+    logging.getLogger("shapesmith.fill").warning("shown")
+    logs.reset()
+    assert "hidden" not in path.read_text() and "shown" in path.read_text()
+
+
+def test_tracebacks_go_to_the_file_only(tmp_path, restore_logging, capsys):
+    path = logs.configure("INFO", tmp_path / "run.log")
+    try:
+        raise RuntimeError("boom")
+    except RuntimeError:
+        logging.getLogger("shapesmith.cli").error("failed: boom", exc_info=True)
+    logs.reset()
+    assert "Traceback" in path.read_text() and "RuntimeError: boom" in path.read_text()
+    console = capsys.readouterr()
+    assert "failed: boom" in console.err and "Traceback" not in console.err and console.out == ""
+
+
+def test_log_path(tmp_path):
+    path = logs.log_path(tmp_path, "ml-export")
+    assert path.parent == tmp_path / "logs" and path.name.startswith("ml-export_") and path.suffix == ".log"
+
+
+def _logging_job(job):
+    logging.getLogger("shapesmith.test_worker").debug(f"worker message {job[0]}")
+    return job[0]
+
+
+def test_worker_messages_reach_the_main_process(tmp_path, restore_logging):
+    path = logs.configure("DEBUG", tmp_path / "run.log")
+    results = sorted(result for _, result, _ in run_jobs(_logging_job, [(1,), (2,)], workers=2, label="test"))
+    logs.reset()
+    text = path.read_text()
+    assert results == [1, 2]
+    assert "worker message 1" in text and "worker message 2" in text
+    assert "MainProcess" not in [line.split()[3] for line in text.splitlines() if "worker message" in line]
+    assert "test: 2/2 jobs done" in text
+
+
+def test_progress_steps():
+    assert logs.progress_steps(3) == [1, 2, 3]
+    assert logs.progress_steps(100) == [10, 20, 30, 40, 50, 60, 70, 80, 90, 100]
+    assert logs.progress_steps(0) == []
+
+
+def test_duration():
+    assert logs.duration(4.24) == "4.2 s"
+    assert logs.duration(125) == "2 min 05 s"
+    assert logs.duration(3725) == "1 h 02 min"