Project

General

Profile

Bug #18025

Need of update translation of GENIE interaction type

Added by Andrea Falcone over 1 year ago. Updated over 1 year ago.

Status:
Work in progress
Priority:
Normal
Category:
Simulation
Target version:
-
Start date:
10/25/2017
Due date:
% Done:

0%

Estimated time:
Spent time:
Occurs In:
Experiment:
ICARUS, SBND
Co-Assignees:
Duration:

Description

While looking at the interaction type of MC neutrino event (generated with icaruscode fcl/gen/genie/simulation_genie_icarus_simple_workshop.fcl), around the 20% of events have an InteractionType()==1000, that is defined as an offset and do not correspond to any physical interaction.
Same error happens also in sbndcode.

GENIEmissingProcessMapping.log.bz2 (114 KB) GENIEmissingProcessMapping.log.bz2 From prodgenie_bnb_nu_uboone.fcl (uboonecode v06_55_01) Gianluca Petrillo, 11/06/2017 02:05 PM
generate_icarus_events.sh (5.16 KB) generate_icarus_events.sh script to generate event similar to those seen to cause issues Robert Hatcher, 11/22/2017 04:17 PM
find_bad_nuance.C (6.92 KB) find_bad_nuance.C script to process events and identify problematic ones and print them out Robert Hatcher, 11/22/2017 04:18 PM
nuance-channels.pdf (638 KB) nuance-channels.pdf nuance channels Robert Hatcher, 11/22/2017 05:11 PM

History

#1 Updated by Gianluca Petrillo over 1 year ago

The module is GENIEGen, the point where this happens is nutools:source:nutools/EventGeneratorBase/GENIE/GENIEHelper.cxx line 1593:

int itype = simb::kNuanceOffset + genie::utils::ghep::NuanceReactionCode(record);

The situation seems to be that nutools stores interaction types in the convention of Nuance, and GENIE provides a helper to convert GENIE codes into Nuance ones.
There may be two cases:
  1. GENIE helper is simulating some interactions for which it does not map the code to Nuance; in this case, GENIE shoudl update their mapping (genie::utils::ghep::NuanceReactionCode())
  2. GENIE is simulating interactions that have no Nuance code; in that case, nutools helper should create specific codes for these GENIE processes

LArSoft code does not get access to the code from GENIE, so we can't inform you on which interaction codes are missing. The FHiCL configuration provided by the reporter should help reproducing the issue, though.

#2 Updated by Robert Hatcher over 1 year ago

Some of GENIE's underlying information is also recorded in a unsullied version in a nusimdata/SimulationBase/GTruth that is associated with every MCTruth object. It would be instructive to see if there is a systematic combination of GTruth.fGint & GTruth.fGscatter that is repeated when this Nuance mapping fails. This would help identify what's going on.

Also, if one could enable GENIEHelper's fcl parameter GHepPrintLevel=13, ala:

https://cdcvs.fnal.gov/redmine/projects/nutools/wiki/GENIEHelper#Printing-GENIE-Event-Record

then we can go back an correlate full GENIE records with MCTruth entries with this nonsense value for simb::MCNeutrino::InteractionType()

For the record, the genie::utils::ghep::NuanceReactionCode(record) code lives here:

https://genie.hepforge.org/trac/browser/generator/trunk/src/GHEP/GHepUtils.cxx#L271

Also a cursory look seems to me that it doesn't properly classify MEC events (but those shouldn't be 20% of a normal sample, I think), so this might be be a "GENIE" issue.

#3 Updated by Andrea Falcone over 1 year ago

I looked for the variable simb::GTruth::fGscatter.
I have found that:
if fGscatter = 10, then InteractionType()=1000,
if InteractionType()=1000 , then fGscatter = 10 or fGscatter = 4.

The event I'm using are in

/pnfs/icarus/scratch/users/cfarnese/workshop/v06_54_00/gen/simulation_genie_icarus_simple_workshop/out/

#4 Updated by Robert Hatcher over 1 year ago

Hmm, that's interesting

Here's GENIE's enumeration of fGscatter

https://genie.hepforge.org/trac/browser/generator/trunk/src/Interaction/ScatteringType.h

4=kScResonant
10=kScMEC,

So, as I suspected GENIE's MEC has no translation to a NUANCE code. This would need a "slot" in MCNeutrino's enumeration to put the value in, and a fix to GENIE's code (and/or a hack to NuTools before GENIE gets a release cut with this fix).

I'm getting a copy of NUANCE and it's scatter type enumeration (soon hopefully) and will look to see where I can start "extending" it to slot in MEC.

The case of fGscatter=4 is a bit more surprising. GENIE's classification of "resonant" scatters are broken down more finely in the NUANCE scheme. GENIE does this translation by looking at the whole event and seeing some topology (#'s of pions, etc). Looks like there might be some case there that isn't handled completely. But to identify the case we would need more full dumps of the event records themselves to figure it out.

This could be done outside of "art"/"larsoft" with basic GENIE functionality. But I'm seriously swamped so if there is someone here at Fermilab who could temporarily be my minion (doing work at my direction) it would happen infinitely faster. They need to be able to setup larsoft (which get's them GENIE), run some gevgen_fnal files (again not hard), and write a bit of C++ as a macro (again not hard).

Volunteers?

#5 Updated by Robert Hatcher over 1 year ago

Also, could I get clarification from the originator which version of GENIE was used, and which set of splines (which genie_xsec UPS product was setup) & the setting of the EventGeneratorList fcl parameter?

For reference:
https://cdcvs.fnal.gov/redmine/projects/nutools/wiki/GENIE_Configuration_Files

#6 Updated by Gianluca Petrillo over 1 year ago

You get the complete set of dependencies with

source /cvmfs/fermilab.opensciencegrid.org/products/larsoft/setup
source /cvmfs/icarus.opensciencegrid.org/products/icarus/setup
ups depend icaruscode v06_54_00 -q e14:prof -f 'Linux64bit+2.6-2.12'
which tells me genie v2_12_8 and genie_xsec v2_12_6. The EventGeneratorList is "Default+CCMEC+NCMEC": the complete GENIEGen module configuration is
      generator: {
         BeamCenter: [
            -1400,
            -350,
            0
         ]
         BeamDirection: [
            0,
            0,
            1
         ]
         BeamName: "booster" 
         BeamRadius: 3
         DebugFlags: 0
         DefinedVtxHistRange: false
         DetectorLocation: "MINOS-NearDet" 
         Environment: [
            "GPRODMODE",
            "YES" 
         ]
         EventGeneratorList: "Default+CCMEC+NCMEC" 
         EventsPerSpill: 1
         FiducialCut: "none" 
         FluxFiles: [
            "converted_beammc_icarus_*.root" 
         ]
         FluxSearchPaths: "/pnfs/icarus/persistent/bnb_gsimple/fluxes_Oct2017/" 
         FluxType: "simple_flux" 
         GHepPrintLevel: 13
         GenFlavors: [
            14
         ]
         GlobalTimeOffset: 10000
         MixerBaseline: 0
         MixerConfig: "none" 
         MonoEnergy: 2
         POTPerSpill: 5e12
         PassEmptySpills: false
         ProductionMode: "yes" 
         RandomTimeOffset: 10000
         SurroundingMass: 0
         TopVolume: "volTPCActive" 
         VtxPosHistRange: [
            0,
            0,
            0,
            0,
            0,
            0
         ]
         XSecTable: "gxspl-FNALsmall.xml" 
         module_type: "GENIEGen" 
      }

#7 Updated by Lynn Garren over 1 year ago

  • Status changed from New to Assigned
  • Assignee set to Robert Hatcher

Gianluca will work with Robert.

#8 Updated by Gianluca Petrillo over 1 year ago

I have instrumented our generator module to dump all the information that LArSoft receives about the generated event, when its interaction type is kNuanceOffset (1000).
Is this information going to be enough? (needless to say, it would make it easier for us to provide the needed information).

I attach the art output from 1000 events generated with a MicroBooNE BNB configuration (GENIEmissingProcessMapping.log.bz2).
For each event, the format is:

%MSG-w GENIEmissingProcessMapping:  GENIEGen:generator@  run: 1 subRun: 0 event: 1
Found an interaction that is not represented by the interaction type code in GENIE:
MCTruth record:
16 particles from neutrinos from beam
neutrino information: from charged weak current, <nuance offset>, mode: 1 (resonant), target: 1000180400 (Ar40), hit nucleon: 2112 (neutron), hit quark: 0
  x=0.136345 y=0.780693 w=1.53629 Q^2=0.233431 GeV^2; theta=0.864647 rad pT=0.00456249 GeV/c
  neutrino:     ID=0  PDG ID 14 mass=0 GeV/c2  status=0 (initial state) weight=0 rescattered at vertex (0.102213, 1.70348, 0.426736; 0)
    created at (282.29, 100.544, 331.583; 4492.84) cm by primary from ID=-1 with momentum (0.00359235, 0.0028127, 1.16783; 1.16784) GeV/c
    still alive!
    comes with a trajectory 0 cm long in 1 points
  outgoing lepton:     ID=4  PDG ID 13 mass=0.10566 GeV/c2  status=1 (stable final state) weight=0 rescattered at vertex (0.102213, 1.70348, 0.426736; 0)
    created at (282.29, 100.544, 331.583; 4492.84) cm by the gods with momentum (0.165356, -0.0653761, 0.151043; 0.256114) GeV/c
    still alive!
    comes with a trajectory 0 cm long in 1 points
[#0] ID=0  PDG ID 14 mass=0 GeV/c2  status=0 (initial state) weight=0 rescattered at vertex (0.102213, 1.70348, 0.426736; 0)
  created at (282.29, 100.544, 331.583; 4492.84) cm by primary from ID=-1 with momentum (0.00359235, 0.0028127, 1.16783; 1.16784) GeV/c
  still alive!
  comes with a trajectory 0 cm long in 1 points:
      (282.29, 100.544, 331.583; 4492.84)
[... all other particles in the record follow]
GENIE truth record:
interaction code: 2, neutrino scattering code: 4 at (2.8229, 1.00544, 3.31583; 4.44156e-08)
probe: nu_mu with cp=(0.124363, -0.216951, 1.21553; 1.24099) hit nucleon with cp=(-0.092137, 0.167659, -0.0363944; 0.929174) GeV (not a sea quark) in target: Ar40 (Z: 18, A: 40)
event interaction weight (genie internal): 1, interaction probability: 2.89285e-13, cross section: 3.22351e-11, differential cross section: 3.9923e-10
particles after reaction, before FSI: 1 pi+, 0 pi-, 1 pi0, 0 p/pbar, 1 n/nbar
  total 3 particles after reaction before FSI: 1/2 charged/neutral, 2 pions, 1 nucleons
process without charmed hadron, resonance: #1
internal (on shell) genie kinematics: Q^2: 0.233431 GeV^2 q^2: -0.233431 GeV^2, w: 1.51909 GeV^2, t: -99999 GeV^2, x: 0.136061, y: 0.760823
FShadSyst: (0, 0, 0; 0)
%MSG

(I have no idea about what FShadSyst might be)

#9 Updated by Robert Hatcher over 1 year ago

Okay, I've done most of what I was hoping a minion would do. But now I need someone to "think" about the results to identify the condition(s) that are failing in the genie::utils::ghep::NanceReactionCode() function's logic.

I've attached to scripts to document and help with the work.

The generate_icarus_events.sh script is how I generated some events. If whomever can't access my /dune/app/users/rhatcher/fix_nuance_codes area and get a copy of icarus.4242.ghep.root, or if they'd like to make more statistics on their own, this is how they could create their own file to work on.

Once a file of events in GENIE format is available then the second script can be invoked via:

   # assumes genie has been setup (ie. dunetpc, icaruscode, larsoft, whatever ...)
   genie -b -q find_bad_nuance.C\(\"icarus.4242.ghep.root\"\)

This script calls the GENIE code for translating event records into Nuance code and if that returns zero, then it tests whether the event is MEC (a known hole in this scheme). If it isn't MEC then it runs a copy of that routine and prints out the GENIE event record so one can look to see how it's failing; i.e. where's the hole in the NuanceReactionCode() logic.

My run of this yielded:

... [lots of printout ...]
Final results: 10000 checked, 1822 MEC, 164 Resonant, 0 other

So, there are 164 Resonant interactions that don't get classified. But I don't have the time to look for the pattern(s) so here I plead for someone else to take this on. The modified copy of the code is in the script and should be fairly easy to read (also to allow you to test whether your "fix" catches everything). if you don't understand how GENIE's event records are formatted let me know and I can help with this.

Once we have a pattern, we can insert code in to nutools to "fix up" this classification until such time that the fix moves upstream to GENIE itself.

Here's an example of a failure:

*** Event #: 9857

|------------------------------------------------------------------------------------------------------------------|
|GENIE GHEP Event Record [print level:   3]                                                                        |
|------------------------------------------------------------------------------------------------------------------|
| Idx |          Name | Ist |        PDG |   Mother  | Daughter  |      Px |      Py |      Pz |       E |      m  |
|------------------------------------------------------------------------------------------------------------------|
|   0 |         nu_mu |   0 |         14 |  -1 |  -1 |   4 |   4 |   0.002 |   0.002 |   1.391 |   1.391 |   0.000 |
|   1 |          Ar40 |   0 | 1000180400 |  -1 |  -1 |   2 |   3 |   0.000 |   0.000 |   0.000 |  37.216 |  37.216 |
|   2 |       neutron |  11 |       2112 |   1 |  -1 |   5 |   5 |  -0.391 |  -0.060 |  -0.179 |   0.927 | **0.940 | M = 0.819
|   3 |          Ar39 |   2 | 1000180390 |   1 |  -1 |  14 |  14 |   0.391 |   0.060 |   0.179 |  36.288 |  36.286 |
|   4 |           mu- |   1 |         13 |   0 |  -1 |  -1 |  -1 |   0.242 |  -0.343 |   0.234 |   0.493 |   0.106 | P = (-0.503,0.714,-0.487)
|   5 |      N+(1440) |   3 |      12212 |   2 |  -1 |   6 |   7 |  -0.630 |   0.285 |   0.978 |   1.826 | **1.440 | M = 1.378
|   6 |        Delta+ |   3 |       2214 |   5 |  -1 |   8 |   9 |  -0.599 |   0.294 |   0.863 |   1.645 |   1.232 |
|   7 |           pi0 |  14 |        111 |   5 |  -1 |  10 |  10 |  -0.032 |  -0.009 |   0.115 |   0.180 |   0.135 | FSI = 3
|   8 |       neutron |  14 |       2112 |   6 |  -1 |  11 |  12 |  -0.456 |  -0.002 |   0.647 |   1.229 |   0.940 | FSI = 4
|   9 |           pi+ |  14 |        211 |   6 |  -1 |  13 |  13 |  -0.142 |   0.295 |   0.216 |   0.417 |   0.140 | FSI = 3
|  10 |           pi0 |   1 |        111 |   7 |  -1 |  -1 |  -1 |   0.094 |   0.010 |   0.076 |   0.181 |   0.135 |
|  11 |       neutron |   1 |       2112 |   8 |  -1 |  -1 |  -1 |  -0.499 |   0.112 |   0.538 |   1.197 |   0.940 |
|  12 |        proton |   1 |       2212 |   8 |  -1 |  -1 |  -1 |   0.125 |  -0.045 |   0.082 |   0.951 |   0.938 |
|  13 |           pi+ |   1 |        211 |   9 |  -1 |  -1 |  -1 |   0.057 |   0.185 |   0.343 |   0.418 |   0.140 |
|  14 |      HadrBlob |  15 | 2000000002 |   3 |  -1 |  -1 |  -1 |  -0.017 |   0.084 |   0.118 |  35.367 | **0.000 | M = 35.366
|------------------------------------------------------------------------------------------------------------------|
|       Fin-Init:                                                |   0.000 |  -0.000 |  -0.000 |  -0.000 |         |
|------------------------------------------------------------------------------------------------------------------|
|       Vertex:          nu_mu @ (x =     1.29272 m, y =     1.43351 m, z =    15.72722 m, t =    1.191694e-07 s)  |
|------------------------------------------------------------------------------------------------------------------|
| Err flag [bits:15->0] : 0000000000000000    |  1st set:                                                     none |
| Err mask [bits:15->0] : 1111111111111111    |  Is unphysical:    NO |   Accepted:   YES                          |
|------------------------------------------------------------------------------------------------------------------|
| sig(Ev) =       7.17431e-39 cm^2  | d2sig(W,Q2;E)/dWdQ2 =     3.35100e-38 cm^2/GeV^3 | Weight =          1.00000 |
|------------------------------------------------------------------------------------------------------------------|

--------------------------------------------------------------------------------------------------------------
GENIE Interaction Summary
--------------------------------------------------------------------------------------------------------------
[-] [Init-State]
 |--> probe        : PDG-code = 14 (nu_mu)
 |--> nucl. target : Z = 18, A = 40, PDG-Code = 1000180400 (Ar40)
 |--> hit nucleon  : PDC-Code = 2112 (neutron)
 |--> hit quark    : no set
 |--> probe 4P     : (E =     1.391272, Px =     0.002316, Py =     0.002085, Pz =     1.391269)
 |--> target 4P    : (E =    37.215526, Px =     0.000000, Py =     0.000000, Pz =     0.000000)
 |--> nucleon 4P   : (E =     0.927105, Px =    -0.390514, Py =    -0.060287, Pz =    -0.178662)
[-] [Process-Info]
 |--> Interaction : Weak[CC]
 |--> Scattering  : RES
[-] [Kinematics]
 |--> *Selected* Bjorken x = 0.365773
 |--> *Selected* Inelasticity y = 0.628192
 |--> *Selected* Momentum transfer Q2 (>0) = 0.707456
 |--> *Selected* Hadronic invariant mass W = 1.377735
[-] [Exclusive Process Info]
 |--> charm prod.  : false |--> strange prod.  : false
 |--> f/s nucleons : N(p) = 0 N(n) = 0
 |--> f/s pions    : N(pi^0) = 0 N(pi^+) = 0 N(pi^-) = 0
 |--> resonance    : P11(1440)
--------------------------------------------------------------------------------------------------------------

FOR CONVENIENCE: IsResonant 1
 np 0 nn 1 npip 1 npim 0 npi0 1
 IsWeakNC() 0 IsWeakCC() 1
 IsNuP() 0 IsNuN() 1
 IsNuBarP() 0 IsNuBarN() 0

#10 Updated by Robert Hatcher over 1 year ago

Also attaching a PDF of the Nuance channel encoding as I got it from Sam Zeller

#11 Updated by Lynn Garren over 1 year ago

Andrea, would you take a look and see if you can identify a pattern in these events?

#12 Updated by Katherine Lato over 1 year ago

  • Status changed from Assigned to Work in progress
  • Assignee set to Andrea Falcone

#13 Updated by Andrea Falcone over 1 year ago

After a first look at the resonant interactions that don't get classified by Robert routine:
- they seem to be "normal" interaction without strange stuff.
- the resonance particle (delta, eta, lambda) decay in pi or gamma -> can the "strange" number of pions in final stat led to the non classification? Strange
I'll look more in details



Also available in: Atom PDF