Skills › Scientific Agent Skills › uncertainty-and-units › format_result.py
format_result.py
327 lines · skills/uncertainty-and-units/scripts/format_result.py at 5f2bc03
1#!/usr/bin/env python3
2"""Round and render a measurement result the way JCGM 100:2008 7.2 requires.
4The uncertainty is rounded to one or two significant digits first, and the
5value is then rounded to that same decimal place. Doing it in the other order,
6or not at all, produces the familiar `12.34567 +/- 0.1` that claims five digits
7of resolution the measurement does not have.
8"""
10from __future__ import annotations
12import argparse
13import math
14import sys
15from decimal import Decimal, InvalidOperation, ROUND_HALF_EVEN, ROUND_HALF_UP
16from typing import Any
18import _common
19from _common import CliError
22ROUNDING_MODES = {"half-even": ROUND_HALF_EVEN, "half-up": ROUND_HALF_UP}
25def _decimal(value: float, *, label: str) -> Decimal:
26 """Convert a float to the shortest exact Decimal that reprs the same."""
28 try:
29 return Decimal(repr(float(value)))
30 except (InvalidOperation, ValueError, OverflowError) as exc:
31 raise CliError(f"{label} is not a finite decimal number") from exc
34def round_to_uncertainty(
35 value: float,
36 uncertainty: float,
37 *,
38 significant_digits: int = 2,
39 rounding: str = "half-even",
40) -> dict[str, Any]:
41 """Round value and uncertainty to a common decimal place.
43 Returns the rounded pair plus the decimal place they share, so callers can
44 render any notation from the same numbers.
45 """
47 if significant_digits not in (1, 2):
48 raise CliError("uncertainty is reported to one or two significant digits")
49 mode = ROUNDING_MODES.get(rounding)
50 if mode is None:
51 raise CliError(f"rounding must be one of: {', '.join(sorted(ROUNDING_MODES))}")
52 if not math.isfinite(value) or not math.isfinite(uncertainty):
53 raise CliError("value and uncertainty must both be finite")
54 if uncertainty < 0:
55 raise CliError("uncertainty must not be negative")
57 exact_value = _decimal(value, label="value")
58 if uncertainty == 0:
59 return {
60 "exact": True,
61 "decimal_place": None,
62 "value": exact_value,
63 "uncertainty": Decimal(0),
64 }
66 exact_uncertainty = _decimal(uncertainty, label="uncertainty")
67 exponent = exact_uncertainty.adjusted()
68 place = exponent - (significant_digits - 1)
69 quantum = Decimal(1).scaleb(place)
70 rounded_uncertainty = exact_uncertainty.quantize(quantum, rounding=mode)
72 # Rounding 0.0996 to two digits gives 0.100, which now carries three
73 # digits. Recompute the place once against the carried value.
74 if rounded_uncertainty.adjusted() != exponent:
75 exponent = rounded_uncertainty.adjusted()
76 place = exponent - (significant_digits - 1)
77 quantum = Decimal(1).scaleb(place)
78 rounded_uncertainty = exact_uncertainty.quantize(quantum, rounding=mode)
80 rounded_value = exact_value.quantize(quantum, rounding=mode)
81 return {
82 "exact": False,
83 "decimal_place": place,
84 "value": rounded_value,
85 "uncertainty": rounded_uncertainty,
86 }
89def _plain(number: Decimal) -> str:
90 """Render a Decimal without exponent notation, keeping trailing zeros."""
92 text = format(number, "f")
93 return "-0" if text == "-0" else text
96# Outside this decade range, positional notation degenerates into a run of
97# zeros, so every rendering switches to scientific form.
98POSITIONAL_RANGE = range(-6, 16)
101def render(rounded: dict[str, Any], unit: str | None) -> dict[str, str]:
102 """Produce the standard notations for an already-rounded pair."""
104 suffix = f" {unit}" if unit else ""
105 latex_unit = f"\\,\\mathrm{{{unit}}}" if unit else ""
106 value = rounded["value"]
107 uncertainty = rounded["uncertainty"]
109 if rounded["exact"]:
110 exponent = value.adjusted() if value != 0 else 0
111 if exponent in POSITIONAL_RANGE:
112 text = f"{_plain(value)}{suffix}"
113 latex = f"${_plain(value)}{latex_unit}$"
114 else:
115 mantissa = _plain(value.scaleb(-exponent))
116 text = f"{mantissa}e{exponent:+03d}{suffix}"
117 latex = (
118 f"${mantissa} \\times 10^{{{exponent}}}{latex_unit}$"
119 )
120 return {
121 "plusminus": f"{text} (zero supplied uncertainty)",
122 "ascii": f"{text} (zero supplied uncertainty)",
123 "parenthetic": f"{text} (zero supplied uncertainty)",
124 "scientific": f"{text} (zero supplied uncertainty)",
125 "latex": latex,
126 }
128 place = rounded["decimal_place"]
129 digits = int(uncertainty.scaleb(-place).to_integral_value())
130 exponent = value.adjusted() if value != 0 else uncertainty.adjusted()
131 mantissa_value = _plain(value.scaleb(-exponent))
132 mantissa_uncertainty = _plain(uncertainty.scaleb(-exponent))
133 scientific = f"({mantissa_value} ± {mantissa_uncertainty})e{exponent:+03d}{suffix}"
134 scientific_concise = f"{mantissa_value}({digits})e{exponent:+03d}{suffix}"
136 # `place > 0` puts the concise digits left of the decimal point, where they
137 # are ambiguous; scientific notation removes the ambiguity.
138 if exponent not in POSITIONAL_RANGE:
139 return {
140 "plusminus": scientific,
141 "ascii": scientific.replace("±", "+/-"),
142 "parenthetic": scientific_concise,
143 "scientific": scientific,
144 "latex": (
145 f"$({mantissa_value} \\pm {mantissa_uncertainty}) "
146 f"\\times 10^{{{exponent}}}{latex_unit}$"
147 ),
148 }
150 value_text = _plain(value)
151 uncertainty_text = _plain(uncertainty)
152 return {
153 "plusminus": f"{value_text} ± {uncertainty_text}{suffix}",
154 "ascii": f"{value_text} +/- {uncertainty_text}{suffix}",
155 "parenthetic": (
156 scientific_concise if place > 0 else f"{value_text}({digits}){suffix}"
157 ),
158 "scientific": scientific,
159 "latex": f"$({value_text} \\pm {uncertainty_text}){latex_unit}$",
160 }
163def build_statement(
164 renderings: dict[str, str],
165 coverage_factor: float | None,
166 coverage_probability: float | None,
167) -> str:
168 """Write the sentence that has to accompany the number."""
170 if coverage_factor is None:
171 return (
172 f"{renderings['plusminus']}, where the stated uncertainty is a combined "
173 "standard uncertainty (coverage factor k = 1)."
174 )
175 probability = ""
176 if coverage_probability is not None:
177 probability = (
178 f", with a caller-supplied coverage probability of approximately "
179 f"{coverage_probability * 100:.0f}%"
180 )
181 return (
182 f"{renderings['plusminus']}, where the stated uncertainty is an expanded "
183 f"uncertainty U = k*u_c with a coverage factor k = {coverage_factor:g}"
184 f"{probability}; the formatter does not verify coverage or the model assumptions."
185 )
188def collect_warnings(
189 value: float, uncertainty: float, significant_digits: int, rounded: dict[str, Any]
190) -> list[str]:
191 """Flag reporting choices that mislead the reader."""
193 warnings: list[str] = []
194 if uncertainty == 0:
195 return ["zero supplied uncertainty does not establish an exact measurement"]
196 leading = _decimal(uncertainty, label="uncertainty").as_tuple().digits[0]
197 if significant_digits == 1 and leading in (1, 2):
198 warnings.append(
199 f"the uncertainty begins with {leading}; rounding it to one significant "
200 "digit changes it by a large fraction, so report two digits"
201 )
202 if value != 0 and uncertainty / abs(value) > 1.0:
203 warnings.append(
204 "the supplied uncertainty exceeds the magnitude of the estimate; "
205 "retain the measured estimate and state the uncertainty type. Whether "
206 "it is consistent with zero depends on the stated interval and model"
207 )
208 if rounded["value"] == 0 and value != 0:
209 warnings.append(
210 "the estimate rounds to zero at the uncertainty's decimal place; the "
211 "measurement does not resolve it from zero"
212 )
213 return warnings
216def build_parser() -> argparse.ArgumentParser:
217 parser = argparse.ArgumentParser(
218 description=(
219 "Round a value and its uncertainty to a common decimal place and render "
220 "the standard notations (JCGM 100:2008 7.2.6)."
221 )
222 )
223 parser.add_argument(
224 "--value", required=True, type=_common.finite_float, help="the estimate"
225 )
226 parser.add_argument(
227 "--uncertainty",
228 required=True,
229 type=_common.non_negative_float,
230 help="standard or expanded uncertainty, in the same unit as the value",
231 )
232 parser.add_argument("--unit", help="unit label appended to every rendering")
233 parser.add_argument(
234 "--significant-digits",
235 type=_common.bounded_int(1, 2),
236 default=2,
237 help="significant digits retained in the uncertainty (default 2)",
238 )
239 parser.add_argument(
240 "--rounding",
241 choices=sorted(ROUNDING_MODES),
242 default="half-even",
243 help="rounding mode (default half-even)",
244 )
245 parser.add_argument(
246 "--coverage-factor",
247 type=_common.positive_float,
248 help="k, when the uncertainty supplied is an expanded uncertainty",
249 )
250 parser.add_argument(
251 "--coverage-probability",
252 type=_common.probability,
253 help="coverage probability associated with k",
254 )
255 parser.add_argument(
256 "--style",
257 choices=("plusminus", "ascii", "parenthetic", "scientific", "latex"),
258 default="plusminus",
259 help="which rendering --format text prints (default plusminus)",
260 )
261 parser.add_argument("--format", choices=("json", "text"), default="json")
262 parser.add_argument("--output", help="write JSON output to this file")
263 parser.add_argument(
264 "--force", action="store_true", help="overwrite an existing output file"
265 )
266 return parser
269def run(arguments: argparse.Namespace) -> dict[str, Any]:
270 """Round, render, and describe one measurement result."""
272 if arguments.coverage_probability is not None and arguments.coverage_factor is None:
273 raise CliError(
274 "a coverage probability only means something alongside --coverage-factor"
275 )
276 rounded = round_to_uncertainty(
277 arguments.value,
278 arguments.uncertainty,
279 significant_digits=arguments.significant_digits,
280 rounding=arguments.rounding,
281 )
282 renderings = render(rounded, arguments.unit)
283 document: dict[str, Any] = {
284 "input": {
285 "value": arguments.value,
286 "uncertainty": arguments.uncertainty,
287 "unit": arguments.unit,
288 "significant_digits": arguments.significant_digits,
289 "rounding": arguments.rounding,
290 },
291 "rounded_value": float(rounded["value"]),
292 "rounded_uncertainty": float(rounded["uncertainty"]),
293 "decimal_place": rounded["decimal_place"],
294 "renderings": renderings,
295 "statement": build_statement(
296 renderings, arguments.coverage_factor, arguments.coverage_probability
297 ),
298 "warnings": collect_warnings(
299 arguments.value,
300 arguments.uncertainty,
301 arguments.significant_digits,
302 rounded,
303 ),
304 }
305 if arguments.value != 0 and arguments.uncertainty > 0:
306 document["relative_uncertainty"] = arguments.uncertainty / abs(arguments.value)
307 return document
310def main(argv: list[str] | None = None) -> int:
311 arguments = build_parser().parse_args(argv)
312 try:
313 document = run(arguments)
314 if arguments.format == "text":
315 print(document["renderings"][arguments.style])
316 else:
317 _common.emit_json(
318 document, output=arguments.output, force=arguments.force
319 )
320 except CliError as error:
321 print(f"error: {error}", file=sys.stderr)
322 return 2
323 return 0
326if __name__ == "__main__":
327 raise SystemExit(main())