nuisance is hosted by Hepforge, IPPP Durham
close Warning: Can't synchronize with repository "(default)" ("(default)" is not readable or not a Git repository.). Look in the Trac log for more information.

Changes between Version 60 and Version 61 of HowToAddSample


Ignore:
Timestamp:
Jan 6, 2017, 5:26:28 PM (7 years ago)
Author:
Clarence Wret
Comment:

--

Legend:

Unmodified
Added
Removed
Modified
  • HowToAddSample

    v60 v61  
    33This is an informal step-by-step guide on how to add samples into the NUISANCE framework. Most of this can be derived by looking at the recently [https://nuisance.hepforge.org/trac/browser/src?order=date&desc=1 implemented measurements] for up-to-date references.
    44
    5 For this tutorial I'll be adding the T2K CC1pi+ H,,2,,O data. The data comes from the [http://t2k-experiment.org/results/nd280data-numu-cc1pi-xs-on-h2o-2015/ T2K site] and [https://arxiv.org/abs/1605.07964 arxiv]. The data is supplied both in a ROOT file and a number of .csv files: I'll be using the ROOT file here.
    6 
    7 After this we'll be able to make data/MC comparisons to the T2K CC1pi+ H,,2,,O data-set for NEUT, GENIE, GiBUU and !NuWro; all in the same framework. We'll get MC predictions, chi2 values, interaction mode contributions, and more using a consistent signal definition, guaranteed the same across the generators.
     5For this tutorial I'll be adding the '''T2K CC1π^+^ H,,2,,O data'''. The data comes from the [http://t2k-experiment.org/results/nd280data-numu-cc1pi-xs-on-h2o-2015/ T2K site] and [https://arxiv.org/abs/1605.07964 arxiv]. The data is supplied both in a ROOT file and a number of .csv files: I'll be using the ROOT file here but very similar steps apply when using .csv or .txt. files.
     6
     7After implementing we'll be able to make data/MC comparisons to the T2K CC1π^+^ H,,2,,O data for NEUT, GENIE, GiBUU and !NuWro, all in the same framework. We'll get MC predictions, χ^2^ values, interaction mode contributions, and more using a consistent signal definition, guaranteed the same across the generators.
    88
    99'''NUISANCE is currently fully written in C++03 and we would appreciate future contributions to adhere to this standard.'''
     
    2424A quick search through the arxiv document points us to reference ![12]. It is also in our [wiki:ExperimentFlux flux list].
    2525
    26 We then generate events in NEUT 5.3.6 with a suitable card-file, see our [wiki:HowToNeut NEUT guide] on how to do this. We need the correct target (water) and flux, and perform the model selections we want (e.g. Rein-Sehgal or Berger-Sehgal coherent model). The procedure is very similar for other generators too.
    27 
    28 Aim to generate around 1M events with all interaction modes turned on. This way we make sure to get all interaction modes into the topologically defined cross-section. We get a small amount of CCQE events which excite a pion from the outgoing nucleon interacting in the nucleus to kick out a pion, leading to a CCQE+1pi ~ CC1pi+ final state, which is signal for this particular sample.
     26I then generate events in NEUT 5.3.6 with a suitable card-file, see our [wiki:HowToNeut NEUT guide] on how to do this. The procedure is very similar for other generators too. '''We need the correct target (water) and flux, and perform the model selections we want (e.g. Rein-Sehgal or Berger-Sehgal coherent model).'''
     27
     28'''Aim to generate around 1M events with all interaction modes turned on.''' This way we make sure to get all interaction modes into the topologically defined cross-section. We get a small amount of CCQE events which excite a pion from the outgoing nucleon interacting in the nucleus to kick out a pion, leading to a CCQE+1π^+^ ~ CC1π^+^ final state, which is signal for this particular sample.
    2929
    3030
    3131== Choosing a cross-section distribution ==
    3232
    33 The T2K CC1pi+ H,,2,,O data release contains various distributions in FIG 4. In this tutorial I'll look at adding one kinematic distribution and one derived distribution: p,,mu,, and E^rec^,,nu,, shown below.
     33The T2K CC1π^+^ H,,2,,O data release contains various distributions in FIG 4. In this tutorial I'll look at adding one kinematic distribution and one derived distribution: p,,µ,, and E^rec^,,ν,, shown below.
    3434
    3535[[Image(T2K_CC1pip_H2O_pmu.png​, 300px)]][[Image(T2K_CC1pip_H2O_enurec.png​, 300px)]]
    3636
    37 In NUISANCE we try to add all available distributions from a publication. '''However''', some distributions will have detector effects "folded" into them (i.e. they will be raw detector-level data). We can only use these if there is some method which bring detector-level variables (e.g. p_mu seen in the detector) to truth-level variable (e.g. p_mu seen after correcting for the detector effects).
     37In NUISANCE we try to add all available distributions from a publication. '''However''', some distributions will have detector effects "folded" into them (i.e. they will be raw detector-level data). We can only use these if there is some method which bring detector-level variables (e.g. p,,µ,, seen in the detector) to truth-level variable (e.g. p,,µ,, seen after correcting for the detector effects).
    3838
    3939
     
    113113== Writing the header ==
    114114
    115 In the case of our T2K CC1pi+ H,,2,,O data, we're dealing with 1D distributions. They should therefore inherit from the `Measurement1D` base class, as mentioned [#point_base earlier]. The `Measurement1D` class is implemented in `src/FitBase/Measurement1D.cxx`.
     115In the case of our T2K CC1π^+^ H,,2,,O data, we're dealing with 1D distributions. They should therefore inherit from the `Measurement1D` base class, as mentioned [#point_base earlier]. The `Measurement1D` class is implemented in `src/FitBase/Measurement1D.cxx`.
    116116
    117117NUISANCE requires a `constructor` and `destructor` for the class, and we'll need to overload `MeasurementBase` methods which define the dependent variable(s) (p,,mu,, in our case) and what our signal is (CC interaction with one muon and one positive pion with no other pions or mesons and any number of nucleons in our case). The `MeasurementBase` functions which we need to overload are `MeasurementBase::isSignal(FitEvent *event)` and `MeasurementBase::FillEventVariables(FitEvent *event)`.
     
    120120
    121121[=#point_signal]
    122 We see that the `src/T2K` directory contains `T2K_SignalDef.cxx`. This extends the `SignalDef` namespace which holds many of the signal definitions for the classes. If you're implementing multiple distributions for the same measurement (e.g. dSig/dp,,mu,, and dSig/dcosmu) we recommend adding the signal definition to the namespace rather than directly in the `isSignal(FitEvent *event)` function in the class implementation. We'll talk about implementing a signal definition [#point_signaldef later].
     122We see that the `src/T2K` directory contains `T2K_SignalDef.cxx`. This extends the `SignalDef` namespace which holds many of the signal definitions for the classes. If you're implementing multiple distributions for the same measurement (e.g. dSig/dp,,µ,, and dSig/dcosθ,,µ,,) we recommend adding the signal definition to the namespace rather than directly in the `isSignal(FitEvent *event)` function in the class implementation. We'll talk about implementing a signal definition [#point_signaldef later].
    123123
    124124
     
    204204If the cross-section is not in cm^2^/nucleon but instead cm^2^/CH (or any other molecule) you will have to scale the `fScaleFactor` by the number of available target nucleons.
    205205
    206 For CCQE, nu,,mu,, interactions only occur off neutrons, which we have 6 of in CH. But the total number of nucleons in CH is 13. So if an experiment reports their cross-section in cm^2^/CH for CCQE nu,,mu,, and our generator output is cm^2^/nucleon we need to do `fScaleFactor*=13./6.`.
    207 
    208 If instead an experiment measures a CC1pi^+^ cross-section in `cm^2^/CH`, we scale by a factor 13 because CC1pi^+^ happens on both neutrons and protons.
     206For CCQE, ν,,µ,, interactions only occur off neutrons, which we have 6 of in CH. But the total number of nucleons in CH is 13. So if an experiment reports their cross-section in cm^2^/CH for CCQE ν,,µ,, and our generator output is cm^2^/nucleon we need to do `fScaleFactor*=13./6.`.
     207
     208If instead an experiment measures a CC1π^+^ cross-section in `cm^2^/CH`, we scale by a factor 13 because CC1π^+^ happens on both neutrons and protons.
    209209}}}
    210210
     
    252252We are now at the point where we can write the `FillEventVariables(FitEvent *event)` implementation mentioned [=#point_header above]. This function defines the dependent variable(s) and how to get them from each `FitEvent` object that gets passed.
    253253
    254 This is where the p,,mu,, and E^rec^,,nu,, implementations differ.
     254This is where the p,,μ,, and E^rec^,,ν,, implementations differ.
    255255
    256256The `FitEvent` class provides numerous getter functions to make `FillEventVariables` as simple as possible.
     
    293293which concludes the `FillEventVariables(FitEvent *event)` function.
    294294
    295 === For the E^rec^,,nu,, distribution ===
    296 
    297 E^rec^,,nu,, needs to be reconstructed from the outgoing particles, so requires a little more work. However it follows very similar steps to the [#point_muon_dist muon case].
    298 
    299 To reconstruct E^rec^,,nu,, we need:
     295=== For the E^rec^,,ν,, distribution ===
     296
     297E^rec^,,ν,, needs to be reconstructed from the outgoing particles, so requires a little more work. However it follows very similar steps to the [#point_muon_dist muon case].
     298
     299To reconstruct E^rec^,,ν,, we need:
    300300
    301301 * Muon momentum
     
    326326}}}
    327327
    328 We should now write a helper function in the `FitUtils` namespace to get `E^rec^,,nu,,` from the muon, pion and neutrino `TLorentzVector`s. [https://arxiv.org/pdf/1605.07964v2.pdf The paper] specifies the exact form in equation 1. The `FitUtils` implementation lives in `src/Utils/FitUtils.cxx`.
     328We should now write a helper function in the `FitUtils` namespace to get `E^rec^,,ν,,` from the muon, pion and neutrino `TLorentzVector`s. [https://arxiv.org/pdf/1605.07964v2.pdf The paper] specifies the exact form in equation 1. The `FitUtils` implementation lives in `src/Utils/FitUtils.cxx`.
    329329
    330330{{{