#include <AssignLabels.h>

Public Member Functions
	AssignLabels ()
	Default constructor. More...

	~AssignLabels ()

InteractionType	GetInteractionType (simb::MCNeutrino &truth)

InteractionType	GetInteractionTypeFromSlice (int nuPDG, bool nuCCNC, int nuMode)

void	GetTopology (const art::Ptr< simb::MCTruth > truth, unsigned int nTopologyHits)

void	PrintTopology ()

unsigned short	GetNProtons ()

unsigned short	GetNPions ()

unsigned short	GetNPizeros ()

unsigned short	GetNNeutrons ()

short	GetPDG ()

unsigned short	TauMode ()

bool	IsAntineutrino ()

unsigned short	GetTopologyType ()

unsigned short	GetTopologyTypeAlt ()

unsigned short	GetProtoDUNEBeamInteractionType (const simb::MCParticle &particle) const

Private Member Functions
unsigned int	GetNeutralDaughterHitsRecursive (const simb::MCParticle &particle) const

int	GetProcessKey (std::string process) const

Private Attributes
unsigned short	nProton

unsigned short	nPion

unsigned short	nPizero

unsigned short	nNeutron

short	pdgCode

unsigned short	tauMode

Detailed Description

Definition at line 20 of file AssignLabels.h.

Constructor & Destructor Documentation

cvn::AssignLabels::AssignLabels ( )

Default constructor.

Definition at line 19 of file AssignLabels.cxx.

   : nProton(0), nPion(0), nPizero(0), nNeutron(0),
     pdgCode(0), tauMode(0)
   {}

cvn::AssignLabels::~AssignLabels ( )

inline

Definition at line 25 of file AssignLabels.h.

25 {};

Member Function Documentation

InteractionType cvn::AssignLabels::GetInteractionType ( simb::MCNeutrino & truth )

Get Interaction_t from pdg, mode and iscc. Setting pdg and mode to zero triggers cosmic ray

Definition at line 26 of file AssignLabels.cxx.

   {
 
     //if(truth->NeutrinoSet())
     //{
     int pdg = truth.Nu().PdgCode();
     bool iscc = truth.CCNC() == simb::kCC;
     int trueMode = truth.Mode();
 
     if(iscc)
     {
       if(abs(pdg) == 14)
       {
         switch(trueMode){
           case simb::kQE: return kNumuQE; break;
           case simb::kRes: return kNumuRes; break;
           case simb::kDIS: return kNumuDIS; break;
           default: return kNumuOther;
         }
       }
       else if(abs(pdg) == 12)
       {
         switch(trueMode){
           case simb::kQE: return kNueQE; break;
           case simb::kRes: return kNueRes; break;
           case simb::kDIS: return kNueDIS; break;
           default: return kNueOther;
         }
       }
       else if(abs(pdg) == 16)
       {
         switch(trueMode){
           case simb::kQE: return kNutauQE; break;
           case simb::kRes: return kNutauRes; break;
           case simb::kDIS: return kNutauDIS; break;
           default: return kNutauOther;
         }
       }
     }
     else if(trueMode==simb::kNuElectronElastic){
       return kNuElectronElastic;
     }
     else return kNC;
     //}
     //else return kCosmic;
 
     return kOther;
   }

InteractionType cvn::AssignLabels::GetInteractionTypeFromSlice	(	int	nuPDG,
		bool	nuCCNC,
		int	nuMode
	)

Definition at line 75 of file AssignLabels.cxx.

   {
 
     if(iscc){
 
       if(abs(pdg) == 14){
         switch(trueMode){
           case simb::kQE: return kNumuQE; break;
           case simb::kRes: return kNumuRes; break;
           case simb::kDIS: return kNumuDIS; break;
           default: return kNumuOther;
         }
       }
       else if(abs(pdg) == 12){
         switch(trueMode){
           case simb::kQE: return kNueQE; break;
           case simb::kRes: return kNueRes; break;
           case simb::kDIS: return kNueDIS; break;
           default: return kNueOther;
         }
       }
       else if(abs(pdg) == 16){
         switch(trueMode){
           case simb::kQE: return kNutauQE; break;
           case simb::kRes: return kNutauRes; break;
           case simb::kDIS: return kNutauDIS; break;
           default: return kNutauOther;
         }
       }
 
     }
     else if(trueMode==simb::kNuElectronElastic){
       return kNuElectronElastic;
     }
     else return kNC;
 
     return kOther;
   }

unsigned int cvn::AssignLabels::GetNeutralDaughterHitsRecursive ( const simb::MCParticle & particle ) const

private

Definition at line 414 of file AssignLabels.cxx.

                                                                                                 {
 
     unsigned int nSimIDEs = 0;
 
     // The backtrack and particle inventory service will be useful here
     art::ServiceHandle<cheat::BackTrackerService> backTrack;
     art::ServiceHandle<cheat::ParticleInventoryService> partService;
 
     for(int d = 0; d < particle.NumberDaughters(); ++d){
 
       const simb::MCParticle *daughter = partService->TrackIdToParticle_P(particle.Daughter(d));     
       unsigned int localSimIDEs = backTrack->TrackIdToSimIDEs_Ps(daughter->TrackId()).size();      
       std::cout << "Got " << localSimIDEs << " hits from " << daughter->PdgCode() << std::endl;
       if(localSimIDEs == 0) localSimIDEs = GetNeutralDaughterHitsRecursive(*daughter);
 
       nSimIDEs += localSimIDEs;
     }
 
     return nSimIDEs;
   }

unsigned short cvn::AssignLabels::GetNNeutrons ( )

inline

Definition at line 37 of file AssignLabels.h.

37 { return nNeutron; };

cvn::AssignLabels::nNeutron

unsigned short nNeutron

Definition: AssignLabels.h:58

unsigned short cvn::AssignLabels::GetNPions ( )

inline

Definition at line 35 of file AssignLabels.h.

35 { return nPion; };

cvn::AssignLabels::nPion

unsigned short nPion

Definition: AssignLabels.h:56

unsigned short cvn::AssignLabels::GetNPizeros ( )

inline

Definition at line 36 of file AssignLabels.h.

36 { return nPizero; };

cvn::AssignLabels::nPizero

unsigned short nPizero

Definition: AssignLabels.h:57

unsigned short cvn::AssignLabels::GetNProtons ( )

inline

Definition at line 34 of file AssignLabels.h.

34 { return nProton; };

cvn::AssignLabels::nProton

unsigned short nProton

Definition: AssignLabels.h:55

short cvn::AssignLabels::GetPDG ( )

inline

Definition at line 38 of file AssignLabels.h.

38 { return pdgCode; };

cvn::AssignLabels::pdgCode

short pdgCode

Definition: AssignLabels.h:59

int cvn::AssignLabels::GetProcessKey ( std::string process ) const

private

Definition at line 371 of file AssignLabels.cxx.

                                                          {
 
   if(process.compare("primary") == 0)                    return 0;
   else if(process.compare("hadElastic") == 0)            return 1;
   else if(process.compare("pi-Inelastic") == 0)          return 2;
   else if(process.compare("pi+Inelastic") == 0)          return 3;
   else if(process.compare("kaon-Inelastic") == 0)        return 4;
   else if(process.compare("kaon+Inelastic") == 0)        return 5;
   else if(process.compare("protonInelastic") == 0)       return 6;
   else if(process.compare("neutronInelastic") == 0)      return 7;
   else if(process.compare("kaon0SInelastic") == 0)       return 8;
   else if(process.compare("kaon0LInelastic") == 0)       return 9;
   else if(process.compare("lambdaInelastic") == 0)       return 10;
   else if(process.compare("omega-Inelastic") == 0)       return 11;
   else if(process.compare("sigma+Inelastic") == 0)       return 12;
   else if(process.compare("sigma-Inelastic") == 0)       return 13;
   else if(process.compare("sigma0Inelastic") == 0)       return 14;
   else if(process.compare("xi-Inelastic") == 0)          return 15;
   else if(process.compare("xi0Inelastic") == 0)          return 16;
   else if(process.compare("anti_protonInelastic") == 0)  return 20;
   else if(process.compare("anti_neutronInelastic") == 0) return 21;
   else if(process.compare("anti_lambdaInelastic") == 0)  return 22;
   else if(process.compare("anti_omega-Inelastic") == 0)  return 23;
   else if(process.compare("anti_sigma+Inelastic") == 0)  return 24;
   else if(process.compare("anti_sigma-Inelastic") == 0)  return 25;
   else if(process.compare("anti_xi-Inelastic") == 0)     return 26;
   else if(process.compare("anti_xi0Inelastic") == 0)     return 27;
 
   else if(process.compare("Decay") == 0)                 return 30;
   else if(process.compare("FastScintillation") == 0)     return 31;
   else if(process.compare("nKiller") == 0)               return 32; // Remove unwanted neutrons: neutron kinetic energy threshold (default 0) or time limit for neutron track
   else if(process.compare("nCapture") == 0)              return 33; // Neutron capture
 
   else if(process.compare("compt") == 0)                 return 40; // Compton Scattering
   else if(process.compare("rayleigh") == 0)              return 41; // Rayleigh Scattering
   else if(process.compare("phot") == 0)                  return 42; // Photoelectric Effect
   else if(process.compare("conv") == 0)                  return 43; // Pair production
   else if(process.compare("CoupledTransportation") == 0) return 44; //
 
   else return -1;
 }

unsigned short cvn::AssignLabels::GetProtoDUNEBeamInteractionType ( const simb::MCParticle & particle ) const

Definition at line 320 of file AssignLabels.cxx.

                                                                                                    {
 
     unsigned short baseProcess = std::numeric_limits<unsigned short>::max();
 
     // The first thing we can do is look at the process key
     std::string processName = particle.EndProcess();
 
     if(GetProcessKey(processName) > -1){
       // Base process gives us a value from 0 to 44
       baseProcess = static_cast<unsigned int>(GetProcessKey(processName));
     }
 
     std::cout << "What interaction type, then? " << processName << std::endl;
 
     // In the case that we have an inelastic interaction, maybe we can do more.
     art::ServiceHandle< cheat::ParticleInventoryService > piService;
 
     unsigned int nPi0 = 0; // Pi-zeros
     unsigned int nPiM = 0; // Pi-minuses
     unsigned int nPiP = 0; // Pi-pluses
     unsigned int nNeu = 0; // Neutrons
     unsigned int nPro = 0; // Protons
     unsigned int nOth = 0; // Everything else
 
     for(int i = 0; i < particle.NumberDaughters(); ++i){
       const simb::MCParticle* daughter = piService->TrackIdToParticle_P(particle.Daughter(i));
       switch(daughter->PdgCode()){
         case 111  : ++nPi0; break;
         case -211 : ++nPiM; break;
         case 211  : ++nPiP; break;
         case 2112 : ++nNeu; break;
         case 2212 : ++nPro; break;
         default   : ++nOth; break;
       }
     }
 
     std::cout << "Base process = " << baseProcess << std::endl;
     std::cout << "Daughters = " << nPi0 << " pi0s, " << nPiM << " pi-s, "
                                 << nPiP << " pi+s, " << nNeu << " neutrons, " 
                                 << nPro << " protons and " << nOth << " other particles." << std::endl;
 
     // If we have a pion with a pi0 in the final state we can flag it as charge exchange
     if(abs(particle.PdgCode()) == 211 && nPi0 == 1){
       return 45; // First free value after those from the truth utility
     }
     else{
       return baseProcess;
     }
   }  

void cvn::AssignLabels::GetTopology	(	const art::Ptr< simb::MCTruth >	truth,
		unsigned int	nTopologyHits = `0`
	)

Definition at line 116 of file AssignLabels.cxx.

                                                                                                {
 
     const simb::MCNeutrino &nu = truth->GetNeutrino();
 
     // Get neutrino flavour
     if(nu.CCNC() == simb::kCC){
       pdgCode = nu.Nu().PdgCode();
     }
     else{
       pdgCode = 1;
     }
 
     // Get tau topology, if necessary
     tauMode = kNotNutau;
     if (abs(pdgCode) == 16) {
       tauMode = kNutauHad;
       for (int p = 0; p < truth->NParticles(); ++p) {
         if (truth->GetParticle(p).StatusCode() != 1) continue;
         int pdg = abs(truth->GetParticle(p).PdgCode());
         int parent = truth->GetParticle(p).Mother();
         while (parent > 0) parent = truth->GetParticle(parent).Mother();
 
         if (parent == 0) {
           if (pdg == 11) {
             tauMode = kNutauE;
             break;
           } else if (pdg == 13) {
             tauMode = kNutauMu;
             break;
           }
         }
       }
     }
 
 //    std::cout << "Neutrino PDG code is " << pdgCode
 //      << ", tau interaction type is " << tauMode << std::endl;
 
     // Now we need to do some final state particle counting.
 //    unsigned int nParticle = truth.NParticles();
 
     nProton = 0;
     nPion = 0; // Charged pions, that is
     nPizero = 0;
     nNeutron = 0;
 
     // We need an instance of the backtracker to find the number of simulated hits for each track
     art::ServiceHandle<cheat::BackTrackerService> backTrack;
     art::ServiceHandle<cheat::ParticleInventoryService> partService;
 
     // Loop over all of the particles
 //    for(unsigned int p = 0; p < nParticle; ++p){
     for(auto const thisPart : partService->MCTruthToParticles_Ps(truth)){
 
 //      const simb::MCParticle& part = truth.GetParticle(p);
       const simb::MCParticle& part = *thisPart;
 
       int pdg = part.PdgCode();
 
       // Make sure this is a final state particle
       if(part.StatusCode() != 1){
         continue;
       }
 
       // Make sure this particle is a daughter of the neutrino
       if(part.Mother() != 0){
         continue;
       }
 
       // GENIE has some fake particles for energy conservation - eg nuclear binding energy. Ignore these
       if(pdg > 2000000000){
         continue;
       }
 
       // Also don't care about nuclear recoils
       if(pdg > 1000000){
         continue;
       }
 
       // Find how many SimIDEs the track has
       unsigned int nSimIDE = backTrack->TrackIdToSimIDEs_Ps(part.TrackId()).size();
 
       // Check if we have more than 100 MeV of kinetic energy
       // float ke = part.E() - part.Mass();
       //    if( ke < 0.0){
       //      continue;
       //    }
 
       // Special case for pi-zeros since it is the decay photons and their pair produced electrons that deposit energy
       if(pdg == 111 || pdg == 2112){
 //        unsigned int nDummyHits = GetNeutralDaughterHitsRecursive(part);
 //        std::cout << "New method of neutral daughters for " << pdg << " = " << nDummyHits << std::endl;
         // Decay photons
         for(int d = 0; d < part.NumberDaughters(); ++d){
           nSimIDE += backTrack->TrackIdToSimIDEs_Ps(part.Daughter(d)).size();
         }
 //        std::cout << "Old method of neutral daughters for " << pdg << " = " << nSimIDE << std::endl;
       }
 
       // Do we pass the number of hits cut?
       if(nSimIDE < nTopologyHits){
         continue;
       }
 
 //      std::cout << "Final state particle " << pdg << " with ke " << ke << " GeV, " << nSimIDE << " true hits and " << part.NumberDaughters() << " daughters" << std::endl;
 
       switch(abs(pdg)){
         case 111 : ++nPizero;  break;
         case 211 : ++nPion;    break;
         case 2112: ++nNeutron; break;
         case 2212: ++nProton;  break;
         default  : break;
       }
 
     }
 
     // Assign the enums based on the counters
     // switch(nProton){
     //   case 0  : top = static_cast<cvn::TopologyType>(top | kTop0proton); break;
     //   case 1  : top = static_cast<cvn::TopologyType>(top | kTop1proton); break;
     //   case 2  : top = static_cast<cvn::TopologyType>(top | kTop2proton); break;
     //   default : top = static_cast<cvn::TopologyType>(top | kTopNproton); break;
     // }
 
     // switch(nPion){
     //   case 0  : top = static_cast<cvn::TopologyType>(top | kTop0pion); break;
     //   case 1  : top = static_cast<cvn::TopologyType>(top | kTop1pion); break;
     //   case 2  : top = static_cast<cvn::TopologyType>(top | kTop2pion); break;
     //   default : top = static_cast<cvn::TopologyType>(top | kTopNpion); break;
     // }
 
     // switch(nPizero){
     //   case 0  : top = static_cast<cvn::TopologyType>(top | kTop0pizero); break;
     //   case 1  : top = static_cast<cvn::TopologyType>(top | kTop1pizero); break;
     //   case 2  : top = static_cast<cvn::TopologyType>(top | kTop2pizero); break;
     //   default : top = static_cast<cvn::TopologyType>(top | kTopNpizero); break;
     // }
 
     // switch(nNeutron){
     //   case 0  : top = static_cast<cvn::TopologyType>(top | kTop0neutron); break;
     //   case 1  : top = static_cast<cvn::TopologyType>(top | kTop1neutron); break;
     //   case 2  : top = static_cast<cvn::TopologyType>(top | kTop2neutron); break;
     //   default : top = static_cast<cvn::TopologyType>(top | kTopNneutron); break;
     // }
 
     std::cout << "Particle counts: " << nProton << ", " << nPion << ", " << nPizero << ", " << nNeutron << std::endl;
 
     // return top;
 
   }

unsigned short cvn::AssignLabels::GetTopologyType ( )

Definition at line 293 of file AssignLabels.cxx.

                                                {
 
     if (abs(pdgCode) == 12) return kTopNue;
     if (abs(pdgCode) == 14) return kTopNumu;
     if (abs(pdgCode) == 16) {
       if (tauMode == kNutauE)   return kTopNutauE;
       if (tauMode == kNutauMu)  return kTopNutauMu;
       if (tauMode == kNutauHad) return kTopNutauHad;
     }
     if (pdgCode == 1) return kTopNC;
     throw std::runtime_error("Topology type not recognised!");
   }

unsigned short cvn::AssignLabels::GetTopologyTypeAlt ( )

Definition at line 306 of file AssignLabels.cxx.

                                                   {
 
     if (abs(pdgCode) == 12) return kTopNueLike;
     if (abs(pdgCode) == 14) return kTopNumuLike;
     if (abs(pdgCode) == 16) {
       if (tauMode == kNutauE)   return kTopNueLike;
       if (tauMode == kNutauMu)  return kTopNumuLike;
       if (tauMode == kNutauHad) return kTopNutauLike;
     }
     if (pdgCode == 1) return kTopNCLike;
     throw std::runtime_error("Topology type not recognised!");
   }

bool cvn::AssignLabels::IsAntineutrino ( )

Definition at line 286 of file AssignLabels.cxx.

                                     {
 
     if (pdgCode < 0) return true;
     else return false;
 
   }

void cvn::AssignLabels::PrintTopology ( )

Definition at line 266 of file AssignLabels.cxx.

                                   {
 
     std::cout << "== Topology Information ==" << std::endl;
 
     std::cout << " - Neutrino PDG code = " << pdgCode << std::endl;
 
     std::cout << " - Number of protons (3 means >2) = " << nProton << std::endl;
 
     std::cout << " - Number of charged pions (3 means >2) = " << nPion << std::endl;
 
     std::cout << " - Number of pizeros (3 means >2) = " << nPizero << std::endl;
 
     std::cout << " - Number of neutrons (3 means >2) = " << nNeutron << std::endl;
 
     std::cout << " - Topology type is " << GetTopologyType() << std::endl;
 
     std::cout << " - Alternate topology type is " << GetTopologyTypeAlt() << std::endl;
 
   }

unsigned short cvn::AssignLabels::TauMode ( )

inline

Definition at line 39 of file AssignLabels.h.

39 { return tauMode; };

cvn::AssignLabels::tauMode

unsigned short tauMode

Definition: AssignLabels.h:60

Member Data Documentation

unsigned short cvn::AssignLabels::nNeutron

private

Definition at line 58 of file AssignLabels.h.

unsigned short cvn::AssignLabels::nPion

private

Definition at line 56 of file AssignLabels.h.

unsigned short cvn::AssignLabels::nPizero

private

Definition at line 57 of file AssignLabels.h.

unsigned short cvn::AssignLabels::nProton

private

Definition at line 55 of file AssignLabels.h.

short cvn::AssignLabels::pdgCode

private

Definition at line 59 of file AssignLabels.h.

unsigned short cvn::AssignLabels::tauMode

private

Definition at line 60 of file AssignLabels.h.

The documentation for this class was generated from the following files:

dunereco/dunereco/CVN/func/AssignLabels.h
dunereco/dunereco/CVN/func/AssignLabels.cxx

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation