Metadata-Version: 2.4
Name: lighthit
Version: 0.2.0a2
Summary: Spectral radiative-transfer response for neutrino-telescope optical modules
Author: Dmitry Naumov
License: BSD-3-Clause
Project-URL: Repository, https://github.com/NeutrinoHit/LightHit
Keywords: radiative transfer,neutrino telescope,Cherenkov,Baikal-GVD
Classifier: Development Status :: 3 - Alpha
Classifier: License :: OSI Approved :: BSD License
Classifier: Intended Audience :: Science/Research
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Topic :: Scientific/Engineering :: Physics
Requires-Python: >=3.11
Description-Content-Type: text/markdown
License-File: LICENSE
Requires-Dist: numpy>=1.26
Requires-Dist: scipy>=1.11
Provides-Extra: accelerate
Requires-Dist: numba>=0.61; extra == "accelerate"
Provides-Extra: notebook
Requires-Dist: jupyterlab>=4; extra == "notebook"
Requires-Dist: ipykernel>=6; extra == "notebook"
Requires-Dist: nbformat>=5; extra == "notebook"
Requires-Dist: nbclient>=0.10; extra == "notebook"
Requires-Dist: nbconvert>=7; extra == "notebook"
Requires-Dist: matplotlib>=3.8; extra == "notebook"
Requires-Dist: plotly>=5; extra == "notebook"
Provides-Extra: viewer
Requires-Dist: plotly>=5; extra == "viewer"
Provides-Extra: dev
Requires-Dist: pytest>=8; extra == "dev"
Requires-Dist: build>=1; extra == "dev"
Dynamic: license-file

# LightHit

Spectral, time-dependent radiative transport and optical-module response.
The production API folds wavelength-dependent water, Cherenkov emission,
measured detector geometry, OM angular acceptance, quantum efficiency and
glass/gel transmission into expected photoelectrons. The medium is homogeneous
and unbounded; scattering uses Henyey–Greenstein with explicit `g`.

## Production API

No global initialization state is used. A `TransportKernel` owns one medium,
one detector and reusable wavelength-dependent caches; all sources and engines
return the same `TransportResponse` contract.

```python
import lighthit as lh

bgvd = lh.load_bgvd_model(
    "/private/path/bgvd-model-master",
    dataset="2021",                 # all 2304 OMs by default
)

kernel = bgvd.kernel(
    lh.KernelConfig(
        cache_directory="/local/cache/lighthit-bgvd",
        wavelength_nodes=9,
        threshold_pe=0.01,
    ),
).build()

laser = lh.IsotropicFlash.monochromatic(
    position_m=[0, 0, 0], photons=1e12, wavelength_nm=532,
)
laser_response = kernel.transport(laser)

track = lh.CherenkovTrack(
    start_m=[0, 0, -50], direction=[0.2, 0.1, 0.97], length_m=120,
)
track_response = kernel.transport(track, method="track")

shower = lh.G4Shower.from_hdf5("event.h5", event=5)
shower_response = kernel.transport(shower, method="axial")
assert shower_response.metadata["detector_angular_model"] == "exact_m_blocks"
```

A method name says how the **source** is handled. How the **module** is handled
is a property of the module: a detector whose angular acceptance has a measured
Legendre bandwidth gets the exact directional kernel of the `m`-block solver,
and the result reports which model was used in
`metadata["detector_angular_model"]`. `method="auto"` chooses `isotropic`,
`track`, or `axial` from the source type; `method="directional"` forces the
exact module model even for an isotropic detector, and the superseded
source-centroid model is kept only as the experimental `axial_centroid` /
`track_centroid` for regression.

`build()` builds the tables that choice will actually need --- call
`kernel.cache_kinds(method)` to see which --- so a directional detector no
longer pays for the multipole tables it will not read. It defaults to the
element-source route; pass `build(method="all")` when both a flash and an event
will be transported. `build(source=...)` resolves both the table *and* the
wavelengths from the source, which matters for a monochromatic flash: it reads
one line, not the quadrature, and building the quadrature for it produced two
tables that were never touched and left the real one to the first transport.

Research variants (`axial_full`, with `generic` as its descriptive alias) share
the same dispatcher but require explicit permission:

```python
control = kernel.transport(
    shower, method="axial_full", allow_experimental=True,
)
print(kernel.available_methods())
```

The BGVD adapter reads the separately supplied private checkout at runtime. It
uses the full `median_om_coordinates_2021.csv`, wavelength-dependent water,
`efficiency(lambda) * transmission_gel_glass(lambda)`, and the private angular
OM polynomial. It deliberately ignores the private scattering indicatrix and
uses HG with `g=0.9`. Neither the checkout nor derived private files are part
of the LightHit distribution.

Cherenkov spectra use two source fields,

```text
S0(lambda) = lambda^-2
S2(lambda) = lambda^-2 * n_phase(lambda)^-2
```

with spatial/angular coefficients stored independently. The complete derivation
and direct mapping from every formula to code is in
[`docs/chapters/12-rte-to-photoelectrons.qmd`](docs/chapters/12-rte-to-photoelectrons.qmd).

Minimal runnable integrations are
`examples/production_laser.py`, `production_track.py`, and
`production_shower.py`.

## Low-level monochromatic solver

`PointGreenSolver` remains the reviewed low-level Green function of an
instantaneous monochromatic point flash and a point isotropic receiver. A
**directed single photon** and an **isotropic single photon** are supported.

## What the calculation returns

`PointGreenSolver.solve` computes the complex response spectrum per unit
effective area, m⁻². The zero-frequency component gives the integrated
signal ("charge"). `result.readout` gives bin integrals, m⁻²; dividing by
the bin width gives the mean registration rate, m⁻²·ns⁻¹.

Detection efficiency is 1, and sensitivity is the same in every direction.
Area is never silently replaced by a bare number: multiplying by a small
effective area gives the expected photon count. There is no OM surface, no
first-entry condition, and no shadowing in this model.

A directed delta flash observed exactly on its own forward ray gives a
singular response; that request is rejected. The API treats
$\cos\theta\ge1-10^{-12}$ as singular. Geometries that close would need a
finite aperture or a finite angular source distribution instead. In the
shipped demonstration the angle is 60°, so the direct light is zero and the
whole signal is scattered. For an isotropic flash the direct light is a
delta function in time with a finite integral; bin integrals across it are
exact.

## Install into a clean checkout

Python 3.11 or newer is required. Commands below run from
`~/Projects/LightHit`:

```bash
python3 --version
python3 -m venv .venv
source .venv/bin/activate
python -m pip install --upgrade pip
python -m pip install -e '.[accelerate,dev]'
python -m pytest -q
```

After publication, the equivalent user installation is:

```bash
python -m pip install 'lighthit[accelerate]'
```

Add the `notebook` or `viewer` extra only when those interfaces are needed.

None of these commands creates a Git commit or publishes anything.

## A first calculation without Jupyter

```bash
python -m lighthit --config examples/point-green.toml \
  --output .build/point-green --plots
```

With the optional acceleration dependency installed, exact first scattering
uses its compiled scalar kernel automatically; select the compiled angular
tail explicitly with `--angular-backend numba`.  Use `--single-backend numpy`
to force the uncompiled reference path for a comparison.

Output:

- `.build/point-green/spectrum.npz`: frequencies and the 0, 1, ≥2 components;
- `.build/point-green/profiles.npz`: unsmeared bins and bins with the given readout;
- `.build/point-green/report.json`: parameters, charge, timings, diagnostics;
- `.build/point-green/figures/`: separate profile and spectrum plots.

Charge without building a time spectrum:

```bash
python -m lighthit --charge-only --repeat 3 --output .build/charge
```

A quick pass and a stricter numerical setting:

```bash
python -m lighthit --preset quick --output .build/quick --plots
python -m lighthit --preset refined --output .build/refined --plots
python scripts/compare_runs.py .build/point-green .build/refined
```

Numerical library threads can be pinned for comparable timings:

```bash
OPENBLAS_NUM_THREADS=1 OMP_NUM_THREADS=1 VECLIB_MAXIMUM_THREADS=1 \
  python -m lighthit --repeat 3 --output .build/benchmark
```

`--repeat` repeats the calculation itself, not loading a saved result. All
timings are kept in the report. This first release has no on-disk transport
cache. Within one `solve` call, the angular solution is shared across every
observation point.

## Notebook

```bash
python -m jupyter lab notebooks/01_point_green.ipynb
```

Jupyter menu: **Run → Run All Cells**. The notebook calls the same package
as the CLI. It shows normalization, geometry, the phase function, the
spectral equations, charge, time profiles, error diagnostics, and timings.

Non-interactive execution and HTML:

```bash
mkdir -p .build/notebook
python -m jupyter nbconvert --to notebook --execute \
  notebooks/01_point_green.ipynb --output 01_point_green.executed \
  --output-dir .build/notebook --ExecutePreprocessor.timeout=600
python -m jupyter nbconvert --to html \
  .build/notebook/01_point_green.executed.ipynb
```

The source notebook is kept in Git without outputs. After running it,
especially against private water, strip outputs before `git add`:

```bash
python scripts/strip_notebook_outputs.py notebooks/01_point_green.ipynb
```

`.gitignore` does not remove outputs already committed to a tracked
notebook.

## Interactive event viewer

The general event viewer shows several G4 showers and point flashes in one
standalone HTML file: 3-D detector geometry, integrated charge by scattering
order, a selectable per-OM time histogram, a full-array time heatmap, and
frame/cumulative animation. Shower axes are drawn as directed spindles and
muons as directed lines with their stored length. Signed Fourier-inversion bins
are displayed with a signed-log colour scale and are never clipped or
renormalized; adjacent bins can be combined ×2/×4 for display.

With the local files under `g4_data`, run all four stored samples (event 5) and
an isotropic laser through the shared fast Numba cache:

```bash
python -m pip install -e '.[accelerate,viewer]'
python scripts/run_event_viewer.py --output .build/event-viewer --threads 4
```

The default triangular array has three clusters. Electron showers are placed
obliquely near individual clusters, the long muons cross the array, and the
laser sits at the triangle centre. Signals are converted to expected
photoelectrons with explicit effective area and efficiency. Before computing a
full spectrum, a cheap axial proxy plus an exact zero-frequency guard skips OMs
below `--threshold-pe 0.01`; set the threshold to zero to disable this.

Open `.build/event-viewer/viewer.html`. The same viewer can load another
portable `viewer.json` through its file button. Input HDF5 files and generated
event arrays remain local and are not package data.

The screening construction and measured speedups are recorded in
`docs/research/signal-screening.md`. The separate cutoff/node/bin-width study
of time ringing is in `docs/research/time-ringing.md`.

## Python API

```python
import numpy as np
from lighthit import PointGreenSolver, SolverSettings, synthetic_medium

medium = synthetic_medium()  # test values close to Baikal water, not a calibration
solver = PointGreenSolver(medium, SolverSettings())
omega = np.linspace(0.0, 1.2, 241)  # rad/ns; include 0 for the charge
result = solver.solve(
    omega,
    displacement_m=[17.32050807568877, 0.0, 10.0],  # detector minus source, m
    direction=[0.0, 0.0, 1.0],
)
print(result.charge_per_m2)
print(result.timings_s)
front = result.front_time_ns[0]
edges = np.arange(front - 30, front + 651, 2.0)
raw = result.readout(edges, sigma_ns=0.0)
measured = result.readout(edges, sigma_ns=3.0)
```

For an isotropic flash, pass `direction=None`. For several detectors, pass
`displacement_m` of shape `(D,3)`. Emission time is set by
`emission_time_ns`; the number of emitted photons by `photons`. A point
isotropic detector needs no normal vector.

## Private BGVD model

The distribution contains **no** private optical tables, OM coefficients,
geometry rows, event files, or derived private responses. `load_bgvd_model`
reads `BaikalWater.py`, `OpticalModule.py`, and the selected geometry CSV from
a separately supplied checkout, recording content hashes rather than copying
the inputs.

```python
import lighthit as lh

bgvd = lh.load_bgvd_model(
    "/full/local/path/to/bgvd-model",
    dataset="2021",
    g=0.9,
)
```

The path may point at a checkout or its `bgvd_model` directory. If the private
package is installed, the path can be omitted. The older
`providers.load_bgvd_water` remains available for a low-level monochromatic
`PointGreenSolver`, but intentionally has no detector response.

Keep private caches and results outside the source tree or under ignored local
directories. The production adapter is checked against a synthetic double and
locally against the separately supplied package; none of that package is
written into the wheel or source distribution.

A run at $g=0.9$ needs its own scan of settings: numbers validated for a
synthetic medium are not automatically a validation of real water, and the
(scattering degree, spatial degree, $k_{\max}$) triple that converges at one
radius does not automatically converge at another. The convergence variables
and their separate roles belong to
[`docs/chapters/15-validation-and-error-budget.qmd`](docs/chapters/15-validation-and-error-budget.qmd).

## What exactly is solved, and where the formulas are

The math is written up as a short Quarto book in `docs/`:

| Chapter | Content |
|---|---|
| `docs/chapters/01-what-lighthit-computes.qmd` | Observable and the complete production route in one chapter |
| `docs/chapters/02-spectral-rte.qmd`–`07-time-and-bins.qmd` | RTE, adjoint Green function, angular solve, spatial and temporal inversion |
| `docs/chapters/08-cherenkov-source.qmd`–`10-directional-om.qmd` | Spectral Cherenkov sources, tracks/showers and directional modules |
| `docs/chapters/11-cache-numba-memory.qmd`–`15-validation-and-error-budget.qmd` | Production cache, end-to-end response, screening, API and validation |
| `docs/appendices/a-notation.qmd`–`f-catalogue.qmd` | Notation, identities, references, rejected methods and reproducibility |

Render it with Quarto (`cd docs && quarto preview`, or `quarto render` for
a static copy in `docs/_book/`); nothing in the book needs private data or
Geant4.

| File | Contents |
|---|---|
| `src/lighthit/medium.py` | Parameters of one spectral node, m and ns |
| `src/lighthit/angular.py` | Tridiagonal adjoint problem, exact free tail |
| `src/lighthit/single.py` | Coordinate-space first order with the full HG function |
| `src/lighthit/single_fast.py` | Optional compiled scalar kernel for the same exact first order |
| `src/lighthit/green.py` | Radial inversion and the combined 0+1+≥2 spectrum |
| `src/lighthit/readout.py` | Bins and instrument smearing, independent of transport |
| `src/lighthit/viewer.py` | Standalone multi-event 3-D viewer, OM charges, time histograms and animation |
| `src/lighthit/providers.py` | Local private-water provider |
| `tests/` | Independent matrix, analytic, and geometric checks |

Every collision order is included. **Orders $\ge2$ use HG coefficients only
up to $L$**, while the exact first order uses the full HG function — an
explicit, stated composite approximation of the full HG response, whose
difference from the full HG answer vanishes as $L\to\infty$. Free angular
transport is not truncated at $L$: its infinite tail is removed
analytically. The spatial angular inversion has its own, independent degree
$J$.

The $k$ range and quadrature, the frequency band, and the frequency step
are all finite as well. `quick`, `balanced`, `refined` change several
spatial/angular settings at once; they do not fix a mathematical error
bound. Unsmeared bins can oscillate right at the light front. The code keeps
negative values and mass before the front; nothing is clipped or renormalized.
The low-level point-source readout integrates orders 0 and 1 directly in
physical time. The production track/shower route bins the ballistic term
directly but Fourier-inverts both scattered components, so both are subject to
the finite frequency grid.
[`docs/chapters/07-time-and-bins.qmd`](docs/chapters/07-time-and-bins.qmd)
keeps these routes separate.

## Model boundaries

The production model includes measured OM orientations, angular acceptance,
spectral efficiency, Cherenkov tracks and showers, wavelength integration and
reusable radial caches. It does not model a finite OM surface, structural
shadowing, boundaries, layered or inhomogeneous media, polarisation, detector
electronics, trigger/noise, or GPU execution.

Sources and attribution for the mathematical construction are recorded in
[PROVENANCE.md](PROVENANCE.md). LightHit is distributed under BSD-3-Clause;
see [LICENSE](LICENSE).
