In binary mode 0 and in reduced binary modes 1-665, AtRIS will generate data for each particle crossing some sensitive detector interface. We designate the interfaces with IIDs (Interface IDs). The interface id 1 corresponds to the interface between the crust (IID=1) and the first air sheet (IID=2). To determine the interface ID of a particle crossing the interface, AtRIS does the following:
At the end of a step, take a look if the beginning of a step was at a boundary
If yes, calculate the dot product of normalized position and momentum vectors
Take the arcus cosine of the resulting value
convert from radians to degrees and round
This last step makes it possible to save the result in form of an 1 byte unsigned char, since the possible values are 0,1,…,180.
In all of the binary modes, we have the following data at our disposal:
Particle type (PDG Code) (4 byte signed int)
Detector Interface ID (DIID) (2 byte unsigned short)
Polar angle in degrees relative to the zenith (1 byte unsigned char)
Particle's kinetic energy at the crossing instance (4 byte single float)
Primary particles vertex energy (4 byte single float)
If we:
./filter simname 13 1
we get all muons crossing the interface between the crust and the first atmospheric sheet. Still this number needs to be normalized:
The simulation was performed with $n_j$ particles in the bin covering the energy range between $E_j$ and $E_{j+1}$.
We can integrate the differential flux to get the real number of particles hitting the surface from within this energy bin: $$ N_j = \int_{\Omega} d\omega \int_S d\sigma\cdot\hat{r} \int^{E_{j+1}}_{E_j} J(E) dE,$$ where:
J, typically expressed in MeV$^{-1}\cdot$cm$^{-3}\cdot$s$^{-1}$, is the differential flux of the primary particles. This integral can be considered separately since almost always we are going to expose a planet to an isotropic particle flux.
$d\sigma$ is a surface element of the detector (the top of the atmosphere)
$S$ is the total area of the detector (the top of the atmosphere)
$d\omega=d\phi d\cos\theta$ is an element of the solid angle
$\hat{R}$ is the unit vector in the direction of $\omega$
$\Omega$, the domain of $\omega$, is a full hemisphere (half a sphere).
For the case of a uniform flux incident from outside on to a sphere, this integral has been solved in the J.D. Soulivan paper:
Since we are shooting particles only from half a hemisphere, $$ \int_{\Omega} d\omega = 2 \pi,$$ while in cases where also particles are allowed to come from the inside, the result would be $4\pi$.
Let us now assume that the incident intensity has an angular dependence $F(\omega)$. Furthermore, in our particular case, this angular dependence will be only due to the polar angle, thus $F=F(\theta)$. When the incident flux is not isotropic, but has such an angular dependence, the proportionality relation between the number of particles (count rate) and the incident intensity is called the gathering power: $$\Gamma_F = \int_\Omega d\omega \int_S d\sigma\cdot\hat{r} F(\omega)= \int_\Omega d\omega F(\omega) \int_S d\sigma\cdot\hat{r}$$.
The directional response function of the instrument (an atmospheric sheet) is determined with: $$A(\omega)=\int_S d\sigma\cdot \hat{r}.$$ With this we get: $$\Gamma_F = \int_\Omega d\omega F(\omega) A(\omega)$$, which can be rewritten as: $$\Gamma_F = 2\pi A \int^{\theta_{MAX}}_{\theta_{MIN}} F(\theta) \cos\theta \sin\theta d\theta$$
If we didn't have the angular dependency of the incident intensity ($F(\theta)$), this would be a perfect example for the application of “integration by substitution”, and the result, which is left for the reader as an exercise, is $1/2 \cdot 2\pi A$, where $A$ is the detector (top of the atmosphere) surface. If the particles are incident from both sides, then the result is $4\pi$.