== - determining the secondary particle flux of a given species ==

=== - Binary data format ===

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)

=== - Normalization ===

If we: <code> ./filter simname 13 1 </code> 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$. 

==== - Analyzing the gathering power of the instrument ====