Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 16 additions & 4 deletions AGENTS.md
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,7 @@ removals. A JOSS paper + citation remain a post-2.0 follow-up.
| `observation.py` | `Observation` / `ObservationTruth` / `Pointing`: the time-series driver, jitter/drift/dither, per-frame truth (phase 1.2). |
| `spectral.py` | Opt-in spectral mode: `QE`, `SED` (relative or absolute via `from_flux_density`), `Spectrum`, `SpectralBandpass`, effective-QE folding, transmission-product helpers (`product`, `from_file`/`from_product`), optional `astropy.units` coercion. |
| `scene/` | The scene/optics layer: `Scene`, `Source` hierarchy (`PointSource`, `ExtendedSource`, `UniformIllumination`, `Catalog`; point/extended sources accept a `flux_sed` absolute SED), PSFs (`GaussianPSF`/`MoffatPSF`/`AiryPSF`/`ArrayPSF`/`EllipticalGaussianPSF`), `Telescope` (+ `Vignetting`/`RadialDistortion`), `Bandpass` (Vega `johnson` / AB `ab` ugriz·Gaia·2MASS + `Extinction`) in `photometry.py`, `Thermal` graybody background in `thermal.py`, `WCSInfo`, `LightCurve`. Renders a photon-rate map; no randomness. |
| `analysis/` | Measurement helpers: `apertures.py` (`aperture_sum`, `centroid`), `ptc.py` (`photon_transfer_curve`). |
| `analysis/` | Measurement helpers: `apertures.py` (`aperture_sum`, `centroid`), `ptc.py` (`photon_transfer_curve`, camera-driven), `characterize.py` (`stack_statistics`, `characterize_dark`, `characterize_flat` — stack-driven, so it runs on *real* detector data as well as simulated frames; `DarkCharacterization.to_config()` returns a `CameraConfig`). |
| `dataset.py` | Scalable raw+truth dataset generation (phase 1.6): `pairs()` → a streaming `PairDataset` (`to_npz`/`to_arrays`), `random_star_fields()` re-iterable scene source. float32-friendly; no global state. |
| `cli.py` | The `getframes` console entry point (phase 1.6): `presets` / `generate` / `dataset` subcommands driven by a TOML config. |
| `presets/` | Preset library. TOML data files in `presets/data/`, loaded via `importlib.resources`. `load_preset`, `available_presets`, `preset_info`. |
Expand All @@ -64,12 +64,24 @@ Data flows one way: `presets` → `CameraConfig` → `Scene` → `Camera` → `b
backend's camera-owned generator. Every generation method accepts a `seed`.
Never call global NumPy/CuPy random state. CPU and GPU streams repeat within
a backend but are statistically, not pixel-for-pixel, matched across devices.
3. **Physics is auditable.** Noise models live as small, documented, pure functions
3. **Fixed patterns come from `fixed_pattern_seed`, never the per-frame RNG.**
Anything that is a property of the *silicon* — PRNU, DSNU, hot pixels, defects,
amplifier gain/offset, structured bias, and the sCMOS per-pixel read-noise RMS —
belongs in `fixed_pattern_maps()` / `FixedPatternMaps`, keyed on
`fixed_pattern_seed`, so it is identical in every frame and therefore removable
by a master frame. Only the *draw* is per-frame (e.g. read noise re-draws the
Gaussian each frame, but its per-pixel sigma does not change). Getting this wrong
is invisible in single-frame spatial statistics and only shows up in per-pixel
statistics *through time* — the failure mode that shipped in
`read_noise_nonuniformity` until it was caught by cross-validation against real
dark stacks. When adding a detector-structure feature, ask whether it should
repeat across frames, and add a test on the temporal statistic if so.
4. **Physics is auditable.** Noise models live as small, documented, pure functions
in `noise.py`. Document the units (electrons vs. ADU) and the model in the
docstring. State assumptions; cite the model form.
4. **Units are explicit.** Field/variable names carry units (`_e`, `_adu`, `_um`,
5. **Units are explicit.** Field/variable names carry units (`_e`, `_adu`, `_um`,
`_c`, `_s`, `_e_per_s`, `_e_per_adu`). Keep this convention.
5. **Typed and validated.** Full type hints (`mypy --strict` passes). Validate
6. **Typed and validated.** Full type hints (`mypy --strict` passes). Validate
inputs in `CameraConfig.__post_init__` and raise informative `ValueError`s.

## Adding a camera preset
Expand Down
72 changes: 72 additions & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,23 @@ to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).

## [Unreleased]

### Fixed

- **sCMOS per-pixel read noise is now a fixed property of the sensor.** The
per-pixel read-noise RMS map implied by `read_noise_nonuniformity` was drawn from
the *per-frame* generator, so it was re-randomised in every frame. Single-frame
spatial statistics were unaffected, but every pixel ended up with the same
expected noise *through time*, which is not how an sCMOS behaves: each pixel has
its own source-follower and column ADC. The map is now built once from
`fixed_pattern_seed` and cached in `FixedPatternMaps`, alongside PRNU and DSNU.
Verified against dark stacks from three real back-illuminated sCMOS cameras
(KURO 1200B, Prime 95B, Marana 4.2B-11): splitting a stack in half and
correlating the two per-pixel temporal-variance maps gives r = 0.89–0.94 on the
real detectors and r = 0.004 with the old model, now r ≈ 0.96.
**This changes generated pixel values** for any configuration with
`read_noise_nonuniformity > 0`; `read_noise_e` and the spatial statistics of a
single frame are unchanged.

### Added

- Full-detector region-of-interest simulation through
Expand All @@ -15,6 +32,61 @@ to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
`CameraConfig.output_resolution`, and active amplifier-boundary properties make
the full-versus-ROI geometry explicit. Exact full-detector split pixels remain
available when an ROI is active.
- **`getframes.analysis.characterize`: detector characterisation from frame
stacks.** Where `photon_transfer_curve` drives a *simulated* camera, this works
on stacks that already exist -- raw data off a real detector, or simulated
frames. `stack_statistics` reduces any iterable of frames (arrays, `Frame`s, a
`dark_series` generator, your own file reader) to per-pixel temporal mean and
variance in one streaming pass, so stacks larger than memory are fine.
`characterize_dark` then measures conversion gain, read noise (with its
per-pixel map, log-normal width and RTS tail), dark current, bias and DSNU from
darks alone -- no flat field needed, because dark charge is Poisson and so
serves as the PTC charge source. `characterize_flat` adds full well, PRNU and
linearity. `DarkCharacterization.to_config()` returns a `CameraConfig`, closing
the loop: measure a real camera, then simulate it. `StackStats.split=True`
additionally gives `temporal_repeatability`, the split-half test that separates
genuine per-pixel noise structure from chi-squared sampling scatter.
New guide (`docs/guides/characterization.md`) and example
(`examples/15_detector_characterization.py`).
- `read_noise_rts_fraction` / `read_noise_rts_factor`: an optional second,
noisier read-noise population modelling the random-telegraph-signal (RTS) pixels
of a real sCMOS array. Measured on three real sensors, ~0.5% of pixels sit above
3x the median read noise where a single log-normal predicts ~0.01%; these are the
pixels that limit faint-source detection. Defaults to off.
- `detector_glow_edge_scale_px`: makes `detector_glow_e_per_s` edge-concentrated
with an exponential falloff, instead of uniform, modelling amplifier glow emitted
at the array periphery. Renormalised so the array mean is unchanged; still fixed
and exposure-scaling, so an exposure-matched master dark removes it. Defaults to
`0` (uniform, the previous behaviour).

### Changed

- The `princeton_instruments_kuro_1200b`, `photometrics_prime_95b`, and
`andor_marana_4_2b_11` presets now carry **measured** conversion gain, read noise,
dark current, bias offset, and non-uniformity terms, fitted from a per-pixel dark
photon-transfer analysis of real frames rather than taken from datasheets. The
largest corrections: conversion gain (1.25-1.3 -> 0.77-0.87 e-/ADU, the low-signal
leg of these dual-gain modes) and `dark_current_nonuniformity` (0.03 -> 0.11-0.33,
which had been roughly an order of magnitude too low). Each preset documents the
operating mode and temperature the values apply to.
- `dark_current_nonuniformity` raised to `0.23` on the remaining sCMOS presets
(`generic_scmos`, `hamamatsu_orca_fusion`, `hamamatsu_orca_quest_2`,
`tucsen_aries_6504_pro`, `andor_cb1_0_5mp`), which previously carried 0.02-0.03 or
omitted the field entirely. `0.23` is the median of the three cameras measured
against real dark stacks (0.11, 0.23, 0.33); each preset documents that it is a
realistic default carried over from characterised hardware rather than a figure
from that camera's datasheet. The same four conventional sCMOS presets also gain
the measured RTS population (`read_noise_rts_fraction = 0.016`, factor 2.65), and
`andor_cb1_0_5mp` / `hamamatsu_orca_quest_2` gain a `read_noise_nonuniformity` of
0.2 where they previously had none at all. `hamamatsu_orca_quest_2` deliberately
keeps no RTS population --- photon-number resolution depends on a tightly screened
read-noise distribution, and importing a tail measured on conventional 11 um sCMOS
would misrepresent it.
- `andor_marana_4_2b_11` gains its measured hot-pixel population
(`hot_pixel_fraction = 1e-4` above 10x the median dark rate).
- `docs/guides/validation.md` documents how to validate a preset against a real dark
stack: measuring conversion gain from darks alone (no flats needed), and the
split-half test for repeatable per-pixel read noise.

## [2.1.1] - 2026-07-26

Expand Down
194 changes: 194 additions & 0 deletions docs/guides/characterization.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,194 @@
# Detector characterisation

`getframes.analysis.characterize` runs the standard bench measurements on stacks
of frames. It takes plain arrays, so it works equally on **frames from a real
detector** and on frames from a simulated [`Camera`][getframes.camera.Camera] —
and the result carries a `to_config()`, so a real camera can be measured, turned
into a [`CameraConfig`][getframes.config.CameraConfig], and then simulated.

```
frames -> stack_statistics per-pixel temporal mean and variance
-> characterize_dark gain, read noise, dark current, bias, DSNU
-> to_config a CameraConfig
-> Camera synthetic frames matching your detector
```

This is the complement to
[`photon_transfer_curve`][getframes.analysis.ptc.photon_transfer_curve], which
drives a *simulated* camera to characterise it. Here the frames come first.

A runnable end-to-end version of everything below is
[`examples/15_detector_characterization.py`](https://github.com/jacotay7/getframes/blob/main/examples/15_detector_characterization.py).

## Step 1 — reduce each stack

Everything is built from one quantity: for each pixel, its mean and variance
*through a stack*. `stack_statistics` computes both in a single streaming pass,
so an iterator over a stack far larger than memory is fine.

```python
from getframes.analysis import stack_statistics

stats = stack_statistics(frames, exposure_s=2.0, split=True)
stats.mean_adu # (h, w) per-pixel temporal mean, ADU
stats.variance_adu2 # (h, w) per-pixel temporal variance, ADU^2
```

`frames` is any iterable of 2-D frames: NumPy arrays, `Frame` objects, a 3-D
cube, a `Camera.dark_series(...)` generator, or your own reader:

```python
def read_raw(path, shape=(1200, 1200)):
"""Stream a flat little-endian uint16 file, one frame at a time."""
n_bytes = shape[0] * shape[1] * 2
with open(path, "rb") as handle:
while chunk := handle.read(n_bytes):
if len(chunk) < n_bytes:
return
yield np.frombuffer(chunk, dtype="<u2").reshape(shape)


stacks = {
exposure: stack_statistics(read_raw(path), exposure_s=exposure)
for exposure, path in my_files.items()
}
```

## Step 2 — characterise from darks

```python
from getframes.analysis import characterize_dark

result = characterize_dark(stacks) # {exposure_s: StackStats}

result.gain_e_per_adu # conversion gain
result.read_noise_e # median per-pixel read noise
result.dark_current_e_per_s # median dark current
result.bias_offset_adu
result.dark_current_nonuniformity # DSNU
result.read_noise_nonuniformity # log-normal width of the read-noise spread
result.read_noise_rts_fraction # pixels above 3x the median (the RTS tail)
result.read_noise_map_e # (h, w) -- the per-pixel maps behind the scalars
```

Use at least three exposures, and make the longest accumulate enough dark charge
to stand clearly above the read noise. Make the shortest as short as the camera
allows: the read noise is measured there.

### Why darks are enough to measure gain

You do not need a flat field. Dark current is itself a Poisson process, so
thermally generated charge is a perfectly good charge source for a photon
transfer curve. For a dark frame,

$$\text{mean}_\text{ADU}(t) = \text{bias} + \frac{Dt}{g}, \qquad
\text{var}_\text{ADU}(t) = \text{RN}_\text{ADU}^2 + \frac{Dt}{g^2}$$

so $\mathrm{d}\,\text{var}/\mathrm{d}\,\text{mean} = 1/g$ and the dark rate $D$
cancels completely.

The mean on its own is degenerate — it only ever tells you $D/g$, and doubling
both leaves every frame identical. What breaks the degeneracy is that Poisson
statistics fix the mean–variance relation *in electrons* with no free parameter,
$\text{var}_e = \text{mean}_e$. Shot noise is an absolute ruler: $\text{SNR} =
\text{mean}/\sqrt{\text{var}} = \sqrt{N}$ is dimensionless and invariant under
rescaling, so it counts discrete charges whatever units you record them in.

`characterize_dark` fits this **per pixel**, which makes it immune to DSNU (each
pixel is its own regression), and fits a *slope* across exposures, which absorbs
the bias pedestal and the read noise into the two intercepts.

!!! warning "Check the Fano factor"
The whole method assumes the dark charge is Poisson. `result.fano_factor`
reports $\text{var}_e/\text{mean}_e$ for the accumulated charge, which should
come out at 1. If it does not, the gain is not trustworthy — suspect a
non-Poisson noise source, a bias step between acquisition sessions, or
saturation.

## Step 3 — rebuild the detector as a config

```python
config = result.to_config(
"my camera",
pixel_size_um=11.0, # things darks cannot see: supply them
full_well_e=80_000.0,
bit_depth=16,
dark_current_ref_temp_c=-20.0, # the temperature the darks were taken at
)
twin = gf.Camera(config)
```

Everything darks can measure is filled in; the rest takes documented placeholders
you should override. `dark_current_ref_temp_c` matters most — stacks carry no
temperature, so without it the config's temperature scaling will be wrong.

## Is your per-pixel noise real, or sampling scatter?

A variance map always looks structured, because estimating a variance from $n$
frames has its own $\chi^2$ scatter. `split=True` gives you the test that tells
them apart: split the stack in half, compute each half's per-pixel variance, and
correlate.

```python
stats = stack_statistics(frames, split=True)
stats.temporal_repeatability # split-half correlation, 0 to 1
stats.fixed_variance_fraction # fraction of the map's spread that is real
```

A detector whose pixels genuinely differ — every sCMOS — gives a high
correlation, because the *same* pixels are noisy in both halves. Uniform noise
gives ~0. Real back-illuminated sCMOS measures **0.89–0.94**.

The most extreme 1% of pixels are excluded before correlating, and on real data
that matters a great deal. A cosmic ray lands in one half only and inflates that
pixel's variance by orders of magnitude, so a handful of them dominate the
covariance: real 60 s Marana darks score **0.006** unclipped against **0.96**
clipped. Use `stats.repeatability(clip_percentile=100.0)` if you want the plain
Pearson correlation.

Note this responds to any fixed per-pixel variance structure, not only read
noise: at long exposures DSNU shows up here too, because a pixel with more dark
current also carries more shot noise.

## Adding flats

Flats measure what darks cannot: full well, PRNU and linearity.

```python
from getframes.analysis import characterize_flat

flat = characterize_flat(flat_stacks, bias_adu=result.bias_offset_adu)
flat.gain_e_per_adu
flat.full_well_e # None if the curve never rolls over
flat.prnu
flat.nonlinearity
```

Sample from near zero up past saturation, and sample densely near the knee.
`full_well_e` comes from the variance peak, which marks the *onset* of
saturation and reads low by roughly the PRNU: the earliest-saturating pixels
start clipping before the array as a whole reaches its ceiling.

Because these are stacks, the variance used is the per-pixel temporal variance,
which is already free of fixed-pattern noise — the usual trick of differencing
flat pairs is unnecessary. PRNU is then measured separately, from the spatial
spread of the time-averaged flat with its shot-noise contribution removed.

## Accuracy

Against a simulated camera with known parameters (72×72, six exposures, 250
frames each — see `tests/test_characterize.py`):

| Parameter | Recovered to |
| --- | --- |
| Conversion gain | 3% |
| Bias offset | 0.1 ADU |
| Dark current | 5% |
| Read noise (median) | 5% |
| DSNU | 15% |
| Read-noise non-uniformity | 15% |
| PRNU (from flats) | 15% |

Accuracy improves with frame count as $1/\sqrt{n}$; the distribution widths need
the most frames, because each pixel's own noise estimate has to be precise
before its spread across pixels is meaningful.
Loading
Loading