Interaction.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*
3  Copyright (c) 2003-2019, The GENIE Collaboration
4  For the full text of the license visit http://copyright.genie-mc.org
5  or see $GENIE/LICENSE
6 
7  Author: Costas Andreopoulos <costas.andreopoulos \at stfc.ac.uk>
8  University of Liverpool & STFC Rutherford Appleton Lab
9 
10  Changes required to implement the GENIE Boosted Dark Matter module
11  were installed by Josh Berger (Univ. of Wisconsin)
12 */
13 //____________________________________________________________________________
14 
15 #include <sstream>
16 
17 #include <TRootIOCtor.h>
18 
25 
26 using namespace genie;
27 using namespace genie::constants;
28 
29 using std::endl;
30 using std::ostringstream;
31 
33 
34 //____________________________________________________________________________
36  ostream & operator<< (ostream& stream, const Interaction & interaction)
37  {
38  interaction.Print(stream);
39  return stream;
40  }
41 }
42 //___________________________________________________________________________
44 TObject()
45 {
46  this->Init();
47 }
48 //___________________________________________________________________________
50 TObject()
51 {
52  this->Init();
53 
54  fInitialState -> Copy (ist);
55  fProcInfo -> Copy (prc);
56 }
57 //___________________________________________________________________________
58 Interaction::Interaction(const Interaction & interaction) :
59 TObject()
60 {
61  this->Init();
62  this->Copy(interaction);
63 }
64 //___________________________________________________________________________
65 Interaction::Interaction(TRootIOCtor*) :
66 TObject(),
67 fInitialState(0),
68 fProcInfo(0),
69 fKinematics(0),
70 fExclusiveTag(0),
71 fKinePhSp(0)
72 {
73 
74 }
75 //___________________________________________________________________________
77 {
78  this->CleanUp();
79 }
80 //___________________________________________________________________________
82 {
83  this->CleanUp();
84  this->Init();
85 }
86 //___________________________________________________________________________
88 {
89  fInitialState = new InitialState ();
90  fProcInfo = new ProcessInfo ();
91  fKinematics = new Kinematics ();
92  fExclusiveTag = new XclsTag ();
93  fKinePhSp = new KPhaseSpace (this);
94 }
95 //___________________________________________________________________________
97 {
98  if ( fInitialState ) delete fInitialState;
99  if ( fProcInfo ) delete fProcInfo;
100  if ( fKinematics ) delete fKinematics;
101  if ( fExclusiveTag ) delete fExclusiveTag;
102  if ( fKinePhSp ) delete fKinePhSp;
103 
104  fInitialState = 0;
105  fProcInfo = 0;
106  fKinematics = 0;
107  fExclusiveTag = 0;
108  fKinePhSp = 0;
109 }
110 //___________________________________________________________________________
111 void Interaction::Copy(const Interaction & interaction)
112 {
113  const InitialState & init = *interaction.fInitialState;
114  const ProcessInfo & proc = *interaction.fProcInfo;
115  const Kinematics & kine = *interaction.fKinematics;
116  const XclsTag & xcls = *interaction.fExclusiveTag;
117 
118  fInitialState -> Copy (init);
119  fProcInfo -> Copy (proc);
120  fKinematics -> Copy (kine);
121  fExclusiveTag -> Copy (xcls);
122 }
123 //___________________________________________________________________________
124 TParticlePDG * Interaction::FSPrimLepton(void) const
125 {
126  int pdgc = this->FSPrimLeptonPdg();
127 
128  if(pdgc) return PDGLibrary::Instance()->Find(pdgc);
129  else return 0;
130 }
131 //___________________________________________________________________________
133 {
134  const ProcessInfo & proc_info = this -> ProcInfo();
135  const InitialState & init_state = this -> InitState();
136 
137  int pdgc = init_state.ProbePdg();
138 
139  LOG("Interaction", pDEBUG) << "Probe PDG code: " << pdgc;
140 
141  if (proc_info.IsNuElectronElastic())
142  return kPdgElectron;
143 
144  // vN (Weak-NC) or eN (EM)
145  if (proc_info.IsWeakNC() || proc_info.IsEM() || proc_info.IsWeakMix() || proc_info.IsDarkMatter()) return pdgc; // EDIT: DM does not change in FS
146 
147  // vN (Weak-CC)
148  else if (proc_info.IsWeakCC()) {
149  int clpdgc;
150  if (proc_info.IsIMDAnnihilation())
151  clpdgc = kPdgMuon;
152  else
153  clpdgc = pdg::Neutrino2ChargedLepton(pdgc);
154  return clpdgc;
155  }
156 
157  LOG("Interaction", pWARN)
158  << "Could not figure out the final state primary lepton pdg code!!";
159 
160  return 0;
161 }
162 //___________________________________________________________________________
163 TParticlePDG * Interaction::RecoilNucleon(void) const
164 {
165  int rnuc = this->RecoilNucleonPdg();
166 
167  if(rnuc) return PDGLibrary::Instance()->Find(rnuc);
168  else return 0;
169 }
170 //___________________________________________________________________________
172 {
173 // Determine the recoil nucleon PDG code
174 
175  const Target & target = fInitialState->Tgt();
176 
177  int recoil_nuc = 0;
178  int struck_nuc = target.HitNucPdg();
179 
181  bool struck_is_nuc = pdg::IsNucleon(struck_nuc);
182  bool is_weak = fProcInfo->IsWeak();
183  bool is_em = fProcInfo->IsEM();
184  bool is_dm = fProcInfo->IsDarkMatter();
185  assert(struck_is_nuc && (is_weak || is_em || is_dm));
186  if(fProcInfo->IsWeakCC()) {
187  recoil_nuc = pdg::SwitchProtonNeutron(struck_nuc); // CC
188  } else {
189  recoil_nuc = struck_nuc; // NC, EM
190  }
191  }
192 
193  if (fProcInfo->IsMEC()) {
194  bool struck_is_2nuc_cluster = pdg::Is2NucleonCluster(struck_nuc);
195  bool is_weak = fProcInfo->IsWeak();
196  bool is_em = fProcInfo->IsEM();
197  assert(struck_is_2nuc_cluster && (is_weak || is_em));
198  if(fProcInfo->IsWeakCC()) {
199  bool isnu = pdg::IsNeutrino(fInitialState->ProbePdg());
200  // nucleon cluster charge should be incremented by +1 for
201  // neutrino CC and by -1 for antineutrino CC
202  int dQ = (isnu) ? +1 : -1;
203  recoil_nuc = pdg::ModifyNucleonCluster(struck_nuc,dQ); // CC
204  }
205  else {
206  recoil_nuc = struck_nuc; // NC, EM
207  }
208  }
209 
210  LOG("Interaction", pDEBUG) << "Recoil nucleon PDG = " << recoil_nuc;
211  return recoil_nuc;
212 }
213 //___________________________________________________________________________
214 void Interaction::SetInitState(const InitialState & init_state)
215 {
217  fInitialState->Copy(init_state);
218 }
219 //___________________________________________________________________________
220 void Interaction::SetProcInfo(const ProcessInfo & proc_info)
221 {
222  if (!fProcInfo) fProcInfo = new ProcessInfo();
223  fProcInfo->Copy(proc_info);
224 }
225 //___________________________________________________________________________
227 {
228  if (!fKinematics) fKinematics = new Kinematics();
229  fKinematics->Copy(kinematics);
230 }
231 //___________________________________________________________________________
232 void Interaction::SetExclTag(const XclsTag & xcls_tag)
233 {
234  if (!fExclusiveTag) fExclusiveTag = new XclsTag();
235  fExclusiveTag->Copy(xcls_tag);
236 }
237 //___________________________________________________________________________
238 string Interaction::AsString(void) const
239 {
240 // Code-ify the interaction in a string to be used as (part of a) cache
241 // branch key.
242 // Template:
243 // nu:x;tgt:x;N:x;q:x(s/v);proc:x;xclv_tag
244 
245  const Target & tgt = fInitialState->Tgt();
246 
247  ostringstream interaction;
248 
249  // If the probe has non-zero mass, then it is DM
250  if (fInitialState->Probe()->PdgCode() == kPdgDarkMatter) {
251  interaction << "dm;";
252  }
253  else {
254  interaction << "nu:" << fInitialState->ProbePdg() << ";";
255  }
256  interaction << "tgt:" << tgt.Pdg() << ";";
257 
258  if(tgt.HitNucIsSet()) {
259  interaction << "N:" << tgt.HitNucPdg() << ";";
260  }
261  if(tgt.HitQrkIsSet()) {
262  interaction << "q:" << tgt.HitQrkPdg()
263  << (tgt.HitSeaQrk() ? "(s)" : "(v)") << ";";
264  }
265 
266  interaction << "proc:" << fProcInfo->InteractionTypeAsString()
267  << "," << fProcInfo->ScatteringTypeAsString() << ";";
268 
269  string xcls = fExclusiveTag->AsString();
270  interaction << xcls;
271  if(xcls.size()>0) interaction << ";";
272 
273  return interaction.str();
274 }
275 //___________________________________________________________________________
276 void Interaction::Print(ostream & stream) const
277 {
278  const string line(110, '-');
279 
280  stream << endl;
281  stream << line << endl;
282 
283  stream << "GENIE Interaction Summary" << endl;
284  stream << line << endl;
285 
286  stream << *fInitialState << endl; // print initial state
287  stream << *fProcInfo; // print process info
288  stream << *fKinematics; // print scattering parameters
289  stream << *fExclusiveTag; // print exclusive process tag
290 
291  stream << line << endl;
292 }
293 //___________________________________________________________________________
295 {
296  this->Copy(interaction);
297  return (*this);
298 }
299 //___________________________________________________________________________
300 //
301 // **** Methods using the "named constructor" C++ idiom ****
302 //
303 //___________________________________________________________________________
305  int target, int probe, ScatteringType_t st, InteractionType_t it)
306 {
307  InitialState init_state (target, probe);
308  ProcessInfo proc_info (st, it);
309 
310  Interaction * interaction = new Interaction(init_state, proc_info);
311  return interaction;
312 }
313 //___________________________________________________________________________
314 Interaction * Interaction::DISCC(int target, int hitnuc, int probe, double E)
315 {
316  Interaction * interaction =
318 
319  InitialState * init_state = interaction->InitStatePtr();
320  init_state->SetProbeE(E);
321  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
322 
323  return interaction;
324 }
325 //___________________________________________________________________________
327  int target, int hitnuc, int hitqrk, bool fromsea, int probe, double E)
328 {
329  Interaction* interaction = Interaction::DISCC(target,hitnuc,probe,E);
330 
331  Target * tgt = interaction->InitStatePtr()->TgtPtr();
332  tgt -> SetHitQrkPdg (hitqrk);
333  tgt -> SetHitSeaQrk (fromsea);
334 
335  return interaction;
336 }
337 //___________________________________________________________________________
339  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
340 {
341  Interaction * interaction =
343 
344  InitialState * init_state = interaction->InitStatePtr();
345  init_state->SetProbeP4(p4probe);
346  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
347 
348  return interaction;
349 }
350 //___________________________________________________________________________
352  int target, int hitnuc, int hitqrk, bool fromsea, int probe,
353  const TLorentzVector & p4probe)
354 {
355  Interaction* interaction = Interaction::DISCC(target,hitnuc,probe,p4probe);
356 
357  Target * tgt = interaction->InitStatePtr()->TgtPtr();
358  tgt -> SetHitQrkPdg (hitqrk);
359  tgt -> SetHitSeaQrk (fromsea);
360 
361  return interaction;
362 }
363 //___________________________________________________________________________
364 Interaction * Interaction::DISNC(int target, int hitnuc, int probe, double E)
365 {
366  Interaction * interaction =
368 
369  InitialState * init_state = interaction->InitStatePtr();
370  init_state->SetProbeE(E);
371  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
372 
373  return interaction;
374 }
375 //___________________________________________________________________________
377  int target, int hitnuc, int hitqrk, bool fromsea, int probe, double E)
378 {
379  Interaction* interaction = Interaction::DISNC(target,hitnuc,probe,E);
380 
381  Target * tgt = interaction->InitStatePtr()->TgtPtr();
382  tgt -> SetHitQrkPdg (hitqrk);
383  tgt -> SetHitSeaQrk (fromsea);
384 
385  return interaction;
386 }
387 //___________________________________________________________________________
389  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
390 {
391  Interaction * interaction =
393 
394  InitialState * init_state = interaction->InitStatePtr();
395  init_state->SetProbeP4(p4probe);
396  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
397 
398  return interaction;
399 }
400 //___________________________________________________________________________
402  int target, int hitnuc, int hitqrk, bool fromsea, int probe,
403  const TLorentzVector & p4probe)
404 {
405  Interaction * interaction = Interaction::DISNC(target,hitnuc,probe,p4probe);
406 
407  Target * tgt = interaction->InitStatePtr()->TgtPtr();
408  tgt -> SetHitQrkPdg (hitqrk);
409  tgt -> SetHitSeaQrk (fromsea);
410 
411  return interaction;
412 }
413 //___________________________________________________________________________
414 Interaction * Interaction::DISEM(int target, int hitnuc, int probe, double E)
415 {
416  Interaction * interaction =
418 
419  InitialState * init_state = interaction->InitStatePtr();
420  init_state->SetProbeE(E);
421  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
422 
423  return interaction;
424 }
425 //___________________________________________________________________________
427  int target, int hitnuc, int hitqrk, bool fromsea, int probe, double E)
428 {
429  Interaction* interaction = Interaction::DISEM(target,hitnuc,probe,E);
430 
431  Target * tgt = interaction->InitStatePtr()->TgtPtr();
432  tgt -> SetHitQrkPdg (hitqrk);
433  tgt -> SetHitSeaQrk (fromsea);
434 
435  return interaction;
436 }
437 //___________________________________________________________________________
439  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
440 {
441  Interaction * interaction =
443 
444  InitialState * init_state = interaction->InitStatePtr();
445  init_state->SetProbeP4(p4probe);
446  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
447 
448  return interaction;
449 }
450 //___________________________________________________________________________
452  int target, int hitnuc, int hitqrk, bool fromsea, int probe,
453  const TLorentzVector & p4probe)
454 {
455  Interaction * interaction = Interaction::DISEM(target,hitnuc,probe,p4probe);
456 
457  Target * tgt = interaction->InitStatePtr()->TgtPtr();
458  tgt -> SetHitQrkPdg (hitqrk);
459  tgt -> SetHitSeaQrk (fromsea);
460 
461  return interaction;
462 }
463 //___________________________________________________________________________
464 Interaction * Interaction::QELCC(int target, int hitnuc, int probe, double E)
465 {
466  Interaction * interaction =
468 
469  InitialState * init_state = interaction->InitStatePtr();
470  init_state->SetProbeE(E);
471  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
472 
473  return interaction;
474 }
475 //___________________________________________________________________________
477  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
478 {
479  Interaction * interaction =
481 
482  InitialState * init_state = interaction->InitStatePtr();
483  init_state->SetProbeP4(p4probe);
484  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
485 
486  return interaction;
487 }
488 //___________________________________________________________________________
489 Interaction * Interaction::QELNC(int target, int hitnuc, int probe, double E)
490 {
491  Interaction * interaction =
493 
494  InitialState * init_state = interaction->InitStatePtr();
495  init_state->SetProbeE(E);
496  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
497 
498  return interaction;
499 }
500 //___________________________________________________________________________
502  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
503 {
504  Interaction * interaction =
506 
507  InitialState * init_state = interaction->InitStatePtr();
508  init_state->SetProbeP4(p4probe);
509  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
510 
511  return interaction;
512 }
513 //___________________________________________________________________________
514 Interaction * Interaction::QELEM(int target, int hitnuc, int probe, double E)
515 {
516  Interaction * interaction =
518 
519  InitialState * init_state = interaction->InitStatePtr();
520  init_state->SetProbeE(E);
521  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
522 
523  return interaction;
524 }
525 //___________________________________________________________________________
527  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
528 {
529  Interaction * interaction =
531 
532  InitialState * init_state = interaction->InitStatePtr();
533  init_state->SetProbeP4(p4probe);
534  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
535 
536  return interaction;
537 }
538 //___________________________________________________________________________
539 Interaction * Interaction::IBD(int target, int hitnuc, int probe, double E)
540 {
541  Interaction * interaction =
543 
544  InitialState * init_state = interaction->InitStatePtr();
545  init_state->SetProbeE(E);
546  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
547 
548  return interaction;
549 }
550 //___________________________________________________________________________
552  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
553 {
554  Interaction * interaction =
556 
557  InitialState * init_state = interaction->InitStatePtr();
558  init_state->SetProbeP4(p4probe);
559  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
560 
561  return interaction;
562 }
563 //___________________________________________________________________________
564 Interaction * Interaction::RESCC(int target, int hitnuc, int probe, double E)
565 {
566  Interaction * interaction =
568 
569  InitialState * init_state = interaction->InitStatePtr();
570  init_state->SetProbeE(E);
571  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
572 
573  return interaction;
574 }
575 //___________________________________________________________________________
577  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
578 {
579  Interaction * interaction =
581 
582  InitialState * init_state = interaction->InitStatePtr();
583  init_state->SetProbeP4(p4probe);
584  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
585 
586  return interaction;
587 }
588 //___________________________________________________________________________
589 Interaction * Interaction::RESNC(int target, int hitnuc, int probe, double E)
590 {
591  Interaction * interaction =
593 
594  InitialState * init_state = interaction->InitStatePtr();
595  init_state->SetProbeE(E);
596  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
597 
598  return interaction;
599 }
600 //___________________________________________________________________________
602  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
603 {
604  Interaction * interaction =
606 
607  InitialState * init_state = interaction->InitStatePtr();
608  init_state->SetProbeP4(p4probe);
609  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
610 
611  return interaction;
612 }
613 //___________________________________________________________________________
614 Interaction * Interaction::RESEM(int target, int hitnuc, int probe, double E)
615 {
616  Interaction * interaction =
617  Interaction::Create(target,probe,kScResonant, kIntEM);
618 
619  InitialState * init_state = interaction->InitStatePtr();
620  init_state->SetProbeE(E);
621  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
622 
623  return interaction;
624 }
625 //___________________________________________________________________________
627  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
628 {
629  Interaction * interaction =
630  Interaction::Create(target,probe,kScResonant, kIntEM);
631 
632  InitialState * init_state = interaction->InitStatePtr();
633  init_state->SetProbeP4(p4probe);
634  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
635 
636  return interaction;
637 }
638 //___________________________________________________________________________
639 Interaction * Interaction::DFRCC(int tgt,int hitnuc, int probe, double E)
640 {
641  Interaction * interaction =
643 
644  InitialState * init_state = interaction->InitStatePtr();
645  init_state->SetProbeE(E);
646  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
647 
648  return interaction;
649 }
650 //___________________________________________________________________________
652  int tgt, int hitnuc, int probe, const TLorentzVector & p4probe)
653 {
654  Interaction * interaction =
656 
657  InitialState * init_state = interaction->InitStatePtr();
658  init_state->SetProbeP4(p4probe);
659  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
660 
661  return interaction;
662 }
663 //___________________________________________________________________________
664 Interaction * Interaction::COHCC(int tgt, int probe, double E)
665 {
666  Interaction * interaction =
668 
669  InitialState * init_state = interaction->InitStatePtr();
670  init_state->SetProbeE(E);
671 
672  return interaction;
673 }
674 //___________________________________________________________________________
676  int tgt, int probe, const TLorentzVector & p4probe)
677 {
678  Interaction * interaction =
680 
681  InitialState * init_state = interaction->InitStatePtr();
682  init_state->SetProbeP4(p4probe);
683 
684  return interaction;
685 }
686 //___________________________________________________________________________
687 Interaction * Interaction::COHNC(int tgt, int probe, double E)
688 {
689  Interaction * interaction =
691 
692  InitialState * init_state = interaction->InitStatePtr();
693  init_state->SetProbeE(E);
694 
695  return interaction;
696 }
697 //___________________________________________________________________________
699  int tgt, int probe, const TLorentzVector & p4probe)
700 {
701  Interaction * interaction =
703 
704  InitialState * init_state = interaction->InitStatePtr();
705  init_state->SetProbeP4(p4probe);
706 
707  return interaction;
708 }
709 //___________________________________________________________________________
710 Interaction * Interaction::COHEl(int tgt, int probe, double E)
711 {
712  Interaction * interaction =
714 
715  InitialState * init_state = interaction->InitStatePtr();
716  init_state->SetProbeE(E);
717 
718  return interaction;
719 }
720 //___________________________________________________________________________
722  int tgt, int probe, const TLorentzVector & p4probe)
723 {
724  Interaction * interaction =
726 
727  InitialState * init_state = interaction->InitStatePtr();
728  init_state->SetProbeP4(p4probe);
729 
730  return interaction;
731 }
732 //___________________________________________________________________________
734 {
735  Interaction * interaction =
737 
738  InitialState * init_state = interaction->InitStatePtr();
739  init_state->SetProbeE(E);
740 
741  return interaction;
742 }
743 //___________________________________________________________________________
744 Interaction * Interaction::IMD(int target, const TLorentzVector & p4probe)
745 {
746  Interaction * interaction =
748 
749  InitialState * init_state = interaction->InitStatePtr();
750  init_state->SetProbeP4(p4probe);
751 
752  return interaction;
753 }
754 //___________________________________________________________________________
755 Interaction * Interaction::AMNuGamma(int tgt, int nuc, int probe, double E)
756 {
757  Interaction * interaction =
759 
760  InitialState * init_state = interaction->InitStatePtr();
761  init_state->SetProbeE(E);
762  init_state->TgtPtr()->SetHitNucPdg(nuc);
763 
764  return interaction;
765 }
766 //___________________________________________________________________________
768  int tgt, int nuc, int probe, const TLorentzVector & p4probe)
769 {
770  Interaction * interaction =
772 
773  InitialState * init_state = interaction->InitStatePtr();
774  init_state->SetProbeP4(p4probe);
775  init_state->TgtPtr()->SetHitNucPdg(nuc);
776 
777  return interaction;
778 }
779 //___________________________________________________________________________
780 Interaction * Interaction::MECCC(int tgt, int ncluster, int probe, double E)
781 {
782  Interaction * interaction =
783  Interaction::Create(tgt, probe, kScMEC, kIntWeakCC);
784 
785  InitialState * init_state = interaction->InitStatePtr();
786  init_state->SetProbeE(E);
787  init_state->TgtPtr()->SetHitNucPdg(ncluster);
788 
789  return interaction;
790 }
791 //___________________________________________________________________________
793  int tgt, int ncluster, int probe, const TLorentzVector & p4probe)
794 {
795  Interaction * interaction =
796  Interaction::Create(tgt, probe, kScMEC, kIntWeakCC);
797 
798  InitialState * init_state = interaction->InitStatePtr();
799  init_state->SetProbeP4(p4probe);
800  init_state->TgtPtr()->SetHitNucPdg(ncluster);
801 
802  return interaction;
803 }
804 //___________________________________________________________________________
805 Interaction * Interaction::MECCC(int tgt, int probe, double E)
806 {
807  Interaction * interaction =
808  Interaction::Create(tgt, probe, kScMEC, kIntWeakCC);
809 
810  InitialState * init_state = interaction->InitStatePtr();
811  init_state->SetProbeE(E);
812 
813  return interaction;
814 }
815 //___________________________________________________________________________
817  int tgt, int probe, const TLorentzVector & p4probe)
818 {
819  Interaction * interaction =
820  Interaction::Create(tgt, probe, kScMEC, kIntWeakCC);
821 
822  InitialState * init_state = interaction->InitStatePtr();
823  init_state->SetProbeP4(p4probe);
824 
825  return interaction;
826 }
827 
828 //___________________________________________________________________________
829 Interaction * Interaction::MECNC(int tgt, int ncluster, int probe, double E)
830 {
831  Interaction * interaction =
832  Interaction::Create(tgt, probe, kScMEC, kIntWeakNC);
833 
834  InitialState * init_state = interaction->InitStatePtr();
835  init_state->SetProbeE(E);
836  init_state->TgtPtr()->SetHitNucPdg(ncluster);
837 
838  return interaction;
839 }
840 //___________________________________________________________________________
842  int tgt, int ncluster, int probe, const TLorentzVector & p4probe)
843 {
844  Interaction * interaction =
845  Interaction::Create(tgt, probe, kScMEC, kIntWeakNC);
846 
847  InitialState * init_state = interaction->InitStatePtr();
848  init_state->SetProbeP4(p4probe);
849  init_state->TgtPtr()->SetHitNucPdg(ncluster);
850 
851  return interaction;
852 }
853 //___________________________________________________________________________
854 Interaction * Interaction::MECEM(int tgt, int ncluster, int probe, double E)
855 {
856  Interaction * interaction =
857  Interaction::Create(tgt, probe, kScMEC, kIntEM);
858 
859  InitialState * init_state = interaction->InitStatePtr();
860  init_state->SetProbeE(E);
861  init_state->TgtPtr()->SetHitNucPdg(ncluster);
862 
863  return interaction;
864 }
865 //___________________________________________________________________________
867  int tgt, int ncluster, int probe, const TLorentzVector & p4probe)
868 {
869  Interaction * interaction =
870  Interaction::Create(tgt, probe, kScMEC, kIntEM);
871 
872  InitialState * init_state = interaction->InitStatePtr();
873  init_state->SetProbeP4(p4probe);
874  init_state->TgtPtr()->SetHitNucPdg(ncluster);
875 
876  return interaction;
877 }
878 //___________________________________________________________________________
879 Interaction * Interaction::GLR(int tgt, double E)
880 {
881  Interaction * interaction =
883 
884  InitialState * init_state = interaction->InitStatePtr();
885  init_state->SetProbeE(E);
886  init_state->TgtPtr()->SetHitNucPdg(0);
887 
888  return interaction;
889 }
890 //___________________________________________________________________________
891 Interaction * Interaction::GLR(int tgt, const TLorentzVector & p4probe)
892 {
893  Interaction * interaction =
895 
896  InitialState * init_state = interaction->InitStatePtr();
897  init_state->SetProbeP4(p4probe);
898  init_state->TgtPtr()->SetHitNucPdg(0);
899 
900  return interaction;
901 }
902 //___________________________________________________________________________
903 Interaction * Interaction::NDecay(int tgt, int decay_mode, int decayed_nucleon)
904 {
905  Interaction * interaction =
907  interaction->ExclTagPtr()->SetDecayMode(decay_mode);
908 
909  InitialState * init_state = interaction->InitStatePtr();
910  init_state->TgtPtr()->SetHitNucPdg(decayed_nucleon);
911 
912  return interaction;
913 }
914 //___________________________________________________________________________
915 Interaction * Interaction::NOsc(int tgt, int annihilation_mode)
916 {
917  Interaction * interaction =
919  interaction->ExclTagPtr()->SetDecayMode(annihilation_mode);
920  return interaction;
921 }
922 //___________________________________________________________________________
923 Interaction * Interaction::ASK(int tgt, int probe, double E)
924 {
925  Interaction * interaction =
927 
928  InitialState * init_state = interaction->InitStatePtr();
929  init_state->SetProbeE(E);
930 
931  return interaction;
932 }
933 //___________________________________________________________________________
935  int tgt, int probe, const TLorentzVector & p4probe)
936 {
937  Interaction * interaction =
939 
940  InitialState * init_state = interaction->InitStatePtr();
941  init_state->SetProbeP4(p4probe);
942 
943  return interaction;
944 }
945 //___________________________________________________________________________
946 Interaction * Interaction::DME(int target, int hitnuc, int probe, double E)
947 {
948  // EDIT: need to be able to create dark matter elastic
949  Interaction * interaction =
951 
952  InitialState * init_state = interaction->InitStatePtr();
953  init_state->SetProbeE(E);
954  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
955 
956  return interaction;
957 }
958 //___________________________________________________________________________
960  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
961 {
962  // EDIT: need to be able to create dark matter elastic
963  Interaction * interaction =
965 
966  InitialState * init_state = interaction->InitStatePtr();
967  init_state->SetProbeP4(p4probe);
968  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
969 
970  return interaction;
971 }
972 //___________________________________________________________________________
973 Interaction * Interaction::DMDI(int target, int hitnuc, int probe, double E)
974 {
975  Interaction * interaction =
977 
978  InitialState * init_state = interaction->InitStatePtr();
979  init_state->SetProbeE(E);
980  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
981 
982  return interaction;
983 }
984 //___________________________________________________________________________
986  int target, int hitnuc, int hitqrk, bool fromsea, int probe, double E)
987 {
988  Interaction* interaction = Interaction::DMDI(target,hitnuc,probe,E);
989 
990  Target * tgt = interaction->InitStatePtr()->TgtPtr();
991  tgt -> SetHitQrkPdg (hitqrk);
992  tgt -> SetHitSeaQrk (fromsea);
993 
994  return interaction;
995 }
996 //___________________________________________________________________________
998  int target, int hitnuc, int probe, const TLorentzVector & p4probe)
999 {
1000  Interaction * interaction =
1002 
1003  InitialState * init_state = interaction->InitStatePtr();
1004  init_state->SetProbeP4(p4probe);
1005  init_state->TgtPtr()->SetHitNucPdg(hitnuc);
1006 
1007  return interaction;
1008 }
1009 //___________________________________________________________________________
1011  int target, int hitnuc, int hitqrk, bool fromsea, int probe,
1012  const TLorentzVector & p4probe)
1013 {
1014  Interaction * interaction = Interaction::DMDI(target,hitnuc,probe,p4probe);
1015 
1016  Target * tgt = interaction->InitStatePtr()->TgtPtr();
1017  tgt -> SetHitQrkPdg (hitqrk);
1018  tgt -> SetHitSeaQrk (fromsea);
1019 
1020  return interaction;
1021 }
static Interaction * COHNC(int tgt, int probe, double E=0)
Basic constants.
bool IsWeak(void) const
static Interaction * IMD(int tgt, double E=0)
bool IsWeakMix(void) const
bool HitSeaQrk(void) const
Definition: Target.cxx:316
bool IsWeakCC(void) const
bool IsNeutrino(int pdgc)
Definition: PDGUtils.cxx:108
const XML_Char * target
Definition: expat.h:268
void SetProbeP4(const TLorentzVector &P4)
set< int >::iterator it
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
static Interaction * DME(int tgt, int nuc, int probe, double E=0)
string ScatteringTypeAsString(void) const
static Interaction * GLR(int tgt, double E=0)
int HitNucPdg(void) const
Definition: Target.cxx:321
TParticlePDG * RecoilNucleon(void) const
recoil nucleon
bool IsNucleon(int pdgc)
Definition: PDGUtils.cxx:309
int RecoilNucleonPdg(void) const
recoil nucleon pdg
int HitQrkPdg(void) const
Definition: Target.cxx:259
InitialState * fInitialState
Initial State info.
Definition: Interaction.h:173
bool IsQuasiElastic(void) const
Definition: ProcessInfo.cxx:67
const int kPdgAntiNuE
Definition: PDGCodes.h:29
const int kPdgDarkMatter
Definition: PDGCodes.h:195
int Pdg(void) const
Definition: Target.h:72
const int kPdgNuMu
Definition: PDGCodes.h:30
Generated/set kinematical variables for an event.
Definition: Kinematics.h:40
int SwitchProtonNeutron(int pdgc)
Definition: PDGUtils.cxx:319
void Copy(const XclsTag &xcls)
copy input XclsTag object
Definition: XclsTag.cxx:148
bool IsInverseBetaDecay(void) const
TParticlePDG * Probe(void) const
static Interaction * IBD(int tgt, int nuc, int probe, double E=0)
static Interaction * NDecay(int tgt, int decay_mode=-1, int decayed_nucleon=0)
void SetProbeE(double E)
bool IsIMDAnnihilation(void) const
void SetKine(const Kinematics &kine)
const int kPdgElectron
Definition: PDGCodes.h:35
XclsTag * fExclusiveTag
Additional info for exclusive channels.
Definition: Interaction.h:176
static Interaction * QELCC(int tgt, int nuc, int probe, double E=0)
static Interaction * MECNC(int tgt, int nuccluster, int probe, double E=0)
string AsString(void) const
static Interaction * RESNC(int tgt, int nuc, int probe, double E=0)
static Interaction * COHCC(int tgt, int probe, double E=0)
Contains minimal information for tagging exclusive processes.
Definition: XclsTag.h:37
int FSPrimLeptonPdg(void) const
final state primary lepton pdg
int ModifyNucleonCluster(int pdgc, int dQ)
Definition: PDGUtils.cxx:327
Summary information for an interaction.
Definition: Interaction.h:56
static Interaction * DMDI(int tgt, int nuc, int probe, double E=0)
static Interaction * RESEM(int tgt, int nuc, int probe, double E=0)
bool IsWeakNC(void) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
bool IsNuElectronElastic(void) const
void Copy(const InitialState &init_state)
ClassImp(Interaction) namespace genie
Definition: Interaction.cxx:32
Float_t E
Definition: plot.C:20
string AsString(void) const
pack into a string code
Definition: XclsTag.cxx:179
static Interaction * MECEM(int tgt, int nuccluster, int probe, double E=0)
static Interaction * QELNC(int tgt, int nuc, int probe, double E=0)
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:44
void Copy(const Kinematics &kine)
Definition: Kinematics.cxx:93
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
Definition: Target.h:41
int ProbePdg(void) const
Definition: InitialState.h:65
Kinematical phase space.
Definition: KPhaseSpace.h:34
void CleanUp(void)
Definition: Interaction.cxx:96
string InteractionTypeAsString(void) const
bool IsDarkMatterElastic(void) const
Definition: ProcessInfo.cxx:72
void Copy(const Interaction &i)
void SetDecayMode(int decay_mode)
Definition: XclsTag.cxx:128
#define pWARN
Definition: Messenger.h:61
TParticlePDG * FSPrimLepton(void) const
final state primary lepton
bool IsMEC(void) const
bool IsEM(void) const
XclsTag * ExclTagPtr(void) const
Definition: Interaction.h:77
bool Is2NucleonCluster(int pdgc)
Definition: PDGUtils.cxx:365
void SetExclTag(const XclsTag &xcls)
Kinematics * fKinematics
kinematical variables
Definition: Interaction.h:175
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:43
static Interaction * ASK(int tgt, int probe, double E=0)
bool HitNucIsSet(void) const
Definition: Target.cxx:300
bool HitQrkIsSet(void) const
Definition: Target.cxx:309
void Print(ostream &stream) const
static Interaction * QELEM(int tgt, int nuc, int probe, double E=0)
enum genie::EScatteringType ScatteringType_t
ostream & operator<<(ostream &stream, const AlgConfigPool &config_pool)
bool IsDarkMatter(void) const
void SetHitNucPdg(int pdgc)
Definition: Target.cxx:188
static Interaction * DISCC(int tgt, int nuc, int probe, double E=0)
Target * TgtPtr(void) const
Definition: InitialState.h:68
TParticlePDG * Find(int pdgc)
Definition: PDGLibrary.cxx:61
static Interaction * DFRCC(int tgt, int nuc, int probe, double E=0)
assert(nhit_max >=nhit_nbins)
int Neutrino2ChargedLepton(int pdgc)
Definition: PDGUtils.cxx:210
InitialState * InitStatePtr(void) const
Definition: Interaction.h:74
const InitialState & InitState(void) const
Definition: Interaction.h:69
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
Interaction & operator=(const Interaction &i)
copy
ProcessInfo * fProcInfo
Process info (scattering, weak current,...)
Definition: Interaction.h:174
const Target & Tgt(void) const
Definition: InitialState.h:67
const int kPdgMuon
Definition: PDGCodes.h:37
static Interaction * RESCC(int tgt, int nuc, int probe, double E=0)
static Interaction * MECCC(int tgt, int nuccluster, int probe, double E=0)
void kinematics()
Definition: kinematics.C:10
static Interaction * COHEl(int tgt, int probe, double E=0)
static Interaction * DISEM(int tgt, int nuc, int probe, double E=0)
static Interaction * Create(int tgt, int probe, ScatteringType_t st, InteractionType_t it)
Most commonly used PDG codes. A set of utility functions to handle PDG codes is provided in PDGUtils...
enum genie::EInteractionType InteractionType_t
void SetInitState(const InitialState &init)
static Interaction * AMNuGamma(int tgt, int nuc, int probe, double E=0)
void SetProcInfo(const ProcessInfo &proc)
void Copy(const ProcessInfo &proc)
KPhaseSpace * fKinePhSp
Kinematic phase space.
Definition: Interaction.h:177
static Interaction * NOsc(int tgt, int annihilation_mode=-1)
static Interaction * DISNC(int tgt, int nuc, int probe, double E=0)
Initial State information.
Definition: InitialState.h:49
#define pDEBUG
Definition: Messenger.h:64