Neutrino interactions with pythia 8310

Never mind the branch name. This uses the latest pythia 8.310.

applications/c8_air_shower.cpp

@@ -147,7 +147,9 @@ int main(int argc, char** argv) {
                    ->needs(opt_Z)
                    ->check(CLI::Range(1, 58))
                    ->group("Primary");
-  app.add_option("-p,--pdg", "PDG code for primary.")
+  app.add_option("-p,--pdg",
+                 "PDG code for primary (p=2212, gamma=22, e-=11, nu_e=12, mu-=13, "
+                 "nu_mu=14, tau=15, nu_tau=16).")
       ->excludes(opt_A)
       ->excludes(opt_Z)
       ->group("Primary");
@@ -397,6 +399,9 @@ int main(int argc, char** argv) {
 
   corsika::pythia8::Decay decayPythia;
 
+  // neutrino interactions
+  corsika::pythia8::NeutrinoInteraction neutrinoPrimaryPythia(E0);
+
   // hadronic photon interactions in resonance region
   corsika::sophia::InteractionModel sophia;
 
@@ -541,9 +546,9 @@ int main(int argc, char** argv) {
   output.add("ZHS", zhs);
 
   // assemble the final process sequence with radio
-  auto sequence =
-      make_sequence(stackInspect, hadronSequence, decayPythia, emCascade, emContinuous,
-                    coreas, zhs, longprof, observationLevel, thinning, cut);
+  auto sequence = make_sequence(stackInspect, neutrinoPrimaryPythia, hadronSequence,
+                                decayPythia, emCascade, emContinuous, coreas, zhs,
+                                longprof, observationLevel, thinning, cut);
 
   /* === END: SETUP PROCESS LIST === */
 
@@ -554,6 +559,7 @@ int main(int argc, char** argv) {
   Cascade EAS(env, tracking, sequence, output, stack);
 
   // print our primary parameters all in one place
+  CORSIKA_LOG_INFO("Primary name: {}", beamCode);
   if (app["--pdg"]->count() > 0) {
     CORSIKA_LOG_INFO("Primary PDG ID:     {}", app["--pdg"]->as<int>());
   } else {
@@ -591,7+597,7 @@
     if (force_interaction) {
       CORSIKA_LOG_INFO("Fixing first interaction at injection point.");
       EAS.forceInteraction();
     }
 
     primaryWriter.recordPrimary(primaryProperties);