Inheritance diagram for NCRadiativeResonant:

Public Member Functions
	NCRadiativeResonant (fhicl::ParameterSet const &p)

	NCRadiativeResonant (NCRadiativeResonant const &)=delete

	NCRadiativeResonant (NCRadiativeResonant &&)=delete

NCRadiativeResonant &	operator= (NCRadiativeResonant const &)=delete

NCRadiativeResonant &	operator= (NCRadiativeResonant &&)=delete

void	cout_stuff (art::Event &e, bool passed)

void	FillTree (art::Event &e, size_t const mct_index, size_t const parent_index, bool const is_nc_delta_radiative)

void	Reset ()

bool	filter (art::Event &e) override

Private Attributes
TTree *	ftree

int	frun

int	fsubrun

int	fevent

int	fnu_pdg

int	fccnc

int	fmode

int	finteraction_type

int	fis_nc_delta_radiative

int	fparent_status_code

int	fparent_pdg

Detailed Description

Definition at line 34 of file NCRadiativeResonant_module.cc.

Constructor & Destructor Documentation

NCRadiativeResonant::NCRadiativeResonant ( fhicl::ParameterSet const & p )

explicit

Definition at line 68 of file NCRadiativeResonant_module.cc.

                                                                     :
   art::EDFilter(p),
     ftree(nullptr) {
 
     art::ServiceHandle<art::TFileService> tfs;
     ftree = tfs->make<TTree>("NCRadiativeResonantFilter", "");
 
     ftree->Branch("run", &frun, "run/I");
     ftree->Branch("subrun", &fsubrun, "subrun/I");
     ftree->Branch("event", &fevent, "event/I");
     ftree->Branch("nu_pdg", &fnu_pdg, "nu_pdg/I");
     ftree->Branch("ccnc", &fccnc, "ccnc/I");
     ftree->Branch("mode", &fmode, "mode/I");
     ftree->Branch("is_nc_delta_radiative", &fis_nc_delta_radiative, "is_nc_delta_radiative/I");
     ftree->Branch("parent_status_code", &fparent_status_code, "parent_status_code/I");
     ftree->Branch("parent_pdg", &fparent_pdg, "parent_pdg/I");
 
 }

NCRadiativeResonant::NCRadiativeResonant ( NCRadiativeResonant const & )

delete

NCRadiativeResonant::NCRadiativeResonant ( NCRadiativeResonant && )

delete

Member Function Documentation

void NCRadiativeResonant::cout_stuff	(	art::Event &	e,
		bool	passed = `false`
	)

Definition at line 89 of file NCRadiativeResonant_module.cc.

                                                                       {
 
   art::ValidHandle<std::vector<simb::MCTruth>> const & ev_mct =
     e.getValidHandle<std::vector<simb::MCTruth>>("generator");
 
   std::cout << passed << "\n"
       << "===========Summary of Truth info===============\n";
   for(simb::MCTruth const & mct : *ev_mct) {
       std::cout<<std::setw(9)<<"track ID ";
     std::cout<<std::setw(12)<<"Pdgcode ";
     std::cout<<std::setw(11)<<"mother ID ";
     std::cout<<std::setw(11)<<"Status code "<<std::endl;
     std::cout << "----------------------------\n";
     for(int i = 0; i < mct.NParticles(); ++i) {
       simb::MCParticle const & mcp = mct.GetParticle(i);
       std::cout<<std::setw(9)<< mcp.TrackId(); 
     std::cout<<std::setw(12)<< mcp.PdgCode();
     std::cout<<std::setw(11)<< mcp.Mother();
     std::cout<<std::setw(11)<< mcp.StatusCode() << "\n";
     }
   }
 
 }

void NCRadiativeResonant::FillTree	(	art::Event &	e,
		size_t const	mct_index,
		size_t const	parent_index,
		bool const	is_nc_delta_radiative
	)

Definition at line 130 of file NCRadiativeResonant_module.cc.

                                           {
 
   Reset();
 
   frun = e.id().run();
   fsubrun = e.id().subRun();
   fevent = e.id().event();
 
   art::ValidHandle<std::vector<simb::MCTruth>> const & ev_mct =
     e.getValidHandle<std::vector<simb::MCTruth>>("generator");
   simb::MCNeutrino const & mcn = ev_mct->at(mct_index).GetNeutrino(); 
 
   fnu_pdg = mcn.Nu().PdgCode();
   fccnc = mcn.CCNC();
   fmode = mcn.Mode();
   finteraction_type = mcn.InteractionType();
 
   fis_nc_delta_radiative = is_nc_delta_radiative;
 
   if(parent_index != SIZE_MAX) {
     fparent_status_code = ev_mct->at(mct_index).GetParticle(parent_index).StatusCode();
     fparent_pdg = ev_mct->at(mct_index).GetParticle(parent_index).PdgCode();
   }    
 
   ftree->Fill();
 
 }

bool NCRadiativeResonant::filter ( art::Event & e )

override

Definition at line 163 of file NCRadiativeResonant_module.cc.

                                              {
 
   art::ValidHandle<std::vector<simb::MCTruth>> const & ev_mct =
     e.getValidHandle<std::vector<simb::MCTruth>>("generator");
 
 
   //looking for resonant photons...
   //  - Particles from NC interaction
   //  - Particles contian photons with StatusCode 1
   //    - the parent of a photon is not delta0 & delta+ ---> pass
   //    - if the parent of a photon is a photon
   //      - the grandparent of a photon is not delta0 & delta+ ---> pass
   for(size_t ith = 0; ith < ev_mct->size(); ++ith) {//loop through MCTruth 
 
     simb::MCTruth const & mct = ev_mct->at(ith);
 //    std::cout<<"Looking at the "<<ith<<"th mctruth particle."<<std::endl;
     if(mct.GetNeutrino().CCNC() != 1) continue;//only look at nc particles;
 
     std::vector<size_t> exiting_photon_parents;
     for(int jth = 0; jth < mct.NParticles(); ++jth) {//loop through MCParticles
       simb::MCParticle const & mcp = mct.GetParticle(jth);
 
       //CHEKC hardcode, TPC filter:
       if(abs(mcp.Vx())>210 ||  abs(mcp.Vy())>210||mcp.Vz()>510 || mcp.Vz()<-1){
       std::cout<<"OUTSIDE TPC x y z ="<<mcp.Vx()<<" "<<mcp.Vy()<<" "<<mcp.Vz()<<std::endl;
       exit(0);
       }
 
       if(mcp.TrackId() != jth) {
         std::cout << "ERROR: " << __LINE__ << " " << __PRETTY_FUNCTION__ << "\nTrackId does not match index\n";
         exit(1);
       }
       if(!(mcp.StatusCode() == 1 && mcp.PdgCode() == 22)) continue;//next, if this is not a photon with StatusCode 1
       exiting_photon_parents.push_back(mcp.Mother());//collect photon's parents
     }
 
     //looking for delta radiative decay, i.e. mother is a delta with statuscode 3 & daughter is a photon. We don't need exiting deltas.
     std::vector<size_t> in_nucleus_photons;
     for(size_t const s : exiting_photon_parents) {
       simb::MCParticle const & mcp = mct.GetParticle(s);
       //2114 delta0; 2214 delta+
       if(abs(mcp.PdgCode()) != 2114 && abs(mcp.PdgCode()) != 2214 ) {//Want mcp as any particles but not delta0+
 //        if(ftree){ 
 //          FillTree(e, ith, mcp.Mother(), true);
 //        }
         cout_stuff(e,true);
         std::cout<<"\nYES a resonant evt: "<<mcp.PdgCode()<<std::endl;
         return true;
       }else if(mcp.PdgCode() == 22) {//collect photons;
         in_nucleus_photons.push_back(mcp.Mother());
       }
     }
 
     for(size_t const s : in_nucleus_photons) {
       simb::MCParticle const & mcp = mct.GetParticle(s);
       if(abs(mcp.PdgCode()) != 2114 && abs(mcp.PdgCode()) != 2214) {
 //        if(ftree) FillTree(e, ith, s, true);
         cout_stuff(e,true);
         std::cout<<"\nYES a resonant evt with 2 photons: "<<mcp.PdgCode()<<std::endl;
         return true;
       }
     }
 
   }
 
 //  if(ftree) FillTree(e, 0, SIZE_MAX, false);
   return false;
 
 }

NCRadiativeResonant& NCRadiativeResonant::operator= ( NCRadiativeResonant const & )

delete

NCRadiativeResonant& NCRadiativeResonant::operator= ( NCRadiativeResonant && )

delete

void NCRadiativeResonant::Reset ( )

Definition at line 114 of file NCRadiativeResonant_module.cc.

                                 {
 
   frun = -1;
   fsubrun = -1;
   fevent = -1;
   fnu_pdg = 0;
   fccnc = -1;
   fmode = -2;
   finteraction_type = -2;
   fis_nc_delta_radiative = -1;
   fparent_status_code = -1; 
   fparent_pdg = 0;
 
 }