CRT: A numerical tool for propagating ultra-high energy cosmic rays through Galactic magnetic field models

Deflection of ultra high energy cosmic rays (UHECRs) by the Galactic magnetic field (GMF) may be sufficiently strong to hinder identification of the UHECR source distribution. A common method for determining the effect of GMF models on source identification efforts is backtracking cosmic rays. We present the public numerical tool CRT for propagating charged particles through Galactic magnetic field models by numerically integrating the relativistic equation of motion. It is capable of both forward- and back-tracking particles with varying compositions through pre-defined and custom user-created magnetic fields. These particles are injected from various types of sources specified and distributed according to the user. Here, we present a description of some source and magnetic field model implementations, as well as validation of the integration routines.

💡 Research Summary

The paper introduces CRT (Cosmic Ray Tracker), an open‑source numerical tool designed to propagate ultra‑high‑energy cosmic rays (UHECRs) through models of the Galactic magnetic field (GMF). The motivation stems from the fact that the GMF can deflect UHECRs by several degrees, severely complicating attempts to associate observed arrival directions with their astrophysical sources. Traditional studies have largely relied on back‑tracking particles from Earth to infer possible source regions, but this approach neglects energy losses, interactions, and the non‑reversible nature of numerical integration in highly structured fields. CRT addresses these shortcomings by supporting both forward and backward propagation, allowing researchers to simulate realistic source distributions, particle compositions, and energy spectra.

At the core of CRT is a relativistic integration of the Lorentz force equation,

 d p/dt = q v × B(r),

implemented with a fourth‑order Runge‑Kutta scheme combined with adaptive step‑size control. This ensures high precision even when particles traverse regions of rapidly varying magnetic field strength. The code dynamically updates particle momentum, charge state, and energy, making it possible to simulate mixed compositions (e.g., protons, helium, carbon, iron) and arbitrary spectral shapes within a single run.

CRT includes several pre‑implemented GMF models that are widely used in the community, such as the Jansson‑Farrar, Sun et al., and Pshirkov models. These models decompose the field into disk, spiral arm, and halo components, each with tunable parameters. In addition, users can supply custom magnetic fields either as analytical expressions or as three‑dimensional grids (FITS, HDF5, etc.). CRT interpolates these grids with trilinear or tricubic schemes, preserving continuity and allowing fine‑scale structures—such as magnetic reversals or localized current sheets—to be represented accurately.

Particle injection is highly flexible. Sources can be defined as point‑like emitters, extended spherical or cylindrical volumes, or via user‑specified probability density functions that control spatial and angular distributions. Energy spectra can follow power‑law, log‑normal, or user‑provided tabulated distributions. Composition fractions are set independently, enabling realistic mixed‑composition simulations that reflect observational hints of a heavy component at the highest energies.

The authors validate the integration routines through two complementary tests. First, they compare CRT trajectories in a uniform magnetic field against analytic circular orbits, finding relative errors below 10⁻⁶ in position and energy conservation. Second, they perform forward and backward propagation through a full Jansson‑Farrar field and demonstrate that the statistical properties of arrival directions, deflection angles, and energy spectra are consistent between the two directions, confirming that numerical irreversibility is negligible for practical purposes.

CRT is written in C++ for performance-critical integration, with a Python wrapper that handles configuration, job submission, and post‑processing. Input parameters are specified in a human‑readable YAML file, and the code supports MPI‑based parallelism, allowing simulations of millions of particles on modern clusters. Output can be written in HDF5, ROOT, or plain CSV formats, and built‑in utilities facilitate visualization with Matplotlib, ParaView, or the yt analysis package.

Application examples illustrate CRT’s scientific utility. Using the Jansson‑Farrar model, the authors back‑track 10⁸ GeV protons and find typical deflections of ~5°, whereas iron nuclei at the same rigidity experience deflections exceeding 20°, highlighting the strong composition dependence of GMF lensing. They also demonstrate forward‑tracking from a synthetic source distribution, showing how the GMF can smear an initially anisotropic sky map into a pattern that mimics isotropy, thereby quantifying the loss of source information.

In summary, CRT provides a comprehensive, validated, and user‑friendly platform for UHECR propagation studies. By offering both forward and backward tracking, support for a wide range of magnetic field configurations, and flexible source/composition specifications, it enables detailed investigations of source identification, magnetic field model testing, and multi‑messenger correlations. The public release, complete documentation, and example scripts lower the barrier for the broader astrophysics community to adopt the tool, and the authors anticipate its integration with upcoming observatories such as AugerPrime, TAx4, and POEMMA to refine our understanding of the Galactic magnetic environment and the origins of the most energetic particles in the universe.