### TABLE OF CONTENTS

- 1. /AnaEvent
- 1.1. AnaEvent/sigma_0pions.dat
- 1.2. AnaEvent/particleIDs
- 1.3. AnaEvent/diff_XXX_dSigma_dTheta_hadron_charge_q.dat
- 1.4. AnaEvent/particleIDs_flag
- 1.5. AnaEvent/set_particleIDs_flag
- 1.6. AnaEvent/event_sigma
- 1.7. AnaEvent/sigma.dat
- 1.8. AnaEvent/event_dSigma_dOmega
- 1.9. AnaEvent/event_dSigma_dE
- 1.10. AnaEvent/event_dSigma_dEcostheta
- 1.11. AnaEvent/event_dSigma_dInvMass
- 1.12. AnaEvent/event_dSigma_dLeptonVariables
- 1.13. AnaEvent/event_dump
- 1.14. AnaEvent/event_INIT
- 1.15. AnaEvent/event_CLEAR
- 1.16. AnaEvent/event_add
- 1.17. AnaEvent/event_getParticle
- 1.18. AnaEvent/makeError_hist
- 1.19. AnaEvent/setError_hist
- 1.20. AnaEvent/event_DumpNumbers
- 1.21. AnaEvent/event_pairPhotons
- 1.22. AnaEvent/event_hadronicEnergy

## /AnaEvent [ Modules ]

NAME

module **AnaEvent**

PURPOSE

This module handles the type "tAnaEvent". It can be used for analyzing your particle vector. See e.g. lowElectronAnalysis.f90 for usage.

This module includes several subroutines to print out cross sections when given an array of events.

The routines given in this module belong to two categories: The first category are the more low-level routines, which act on building up an single tAnaEvent or dumping a list of tAnaEvent types, as e.g.:

- event_init
- event_clear
- event_add
- event_dump

On the other hand we have routines as:

which fill histograms by looping over a given list of tAnaEvent types. These are very complex.

INPUTS

No Namelist available.

NOTES

%perweight should not include any dOmega or the such. The sum of all perweight should simply return the total cross section.

## AnaEvent/diff_XXX_dSigma_dTheta_hadron_charge_q.dat [ Output files ]

[ Top ] [ AnaEvent ] [ Output files ]

NAME

file **diff_XXX_dSigma_dTheta_hadron_charge_q.dat**

PURPOSE

The file is produced in the runs with eventtype=5=neutrino .

The file shows the angular distributions after final state interactions for various mesons and baryons with charge q The cross sections are given for events with 1 meson of a given kind, but there can be other mesons or baryons around

Units:

- For process_ID=CC and NC: 10^{-38} cm^2/rad
- For process_ID=EM: nanobarns=10^{-33}cm^2/ rad
- All cross sections are given per nucleon (1/A)

Columns:

- #1: angle in rad (runs from 0 to pi)
- #2: dsigma/dtheta in 10^{-38} cm^2/rad
- #3: number of events in that bin
- #4: statistical error

HERE XXX=000,QE, Delta, denotes the subprocesses, with 000 being the sum of all

## AnaEvent/particleIDs_flag [ Global module-variables ]

[ Top ] [ AnaEvent ] [ Global module-variables ]

SOURCE

` logical, dimension(1:numStableParts) :: `**particleIDs_flag** = .true.

PURPOSE

Flags to switch on/off output of special ID's,
e.g. **particleIDs_flag**(1)=.false. means that there is no output for pions
The array particleIDs is set in module AnaEventDefinition.f90

## AnaEvent/set_particleIDs_flag [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **set_particleIDs_flag**(flags)

PURPOSE

Set the variable **set_particleIDs_flag** to a new value

INPUTS

- logical, dimension(1:numStableParts) :: flags

OUTPUT

- sets particleIds_flag=flags

## AnaEvent/event_sigma [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_sigma**(E,sigma,newInit,runNumber,identifier)

PURPOSE

- Output for single particle production sigma(total), e.g. pi+, pi-
- Output for double particle production sigma(total), e.g. 2pi, 2N
- Output for three particle production sigma(total), e.g. 3pi, 3N
- Output for four particle production sigma(total), e.g. 4pi, 4N

INPUTS

- type(tAnaEvent), dimension(:) :: E -- List of events
- real, dimension(1:dimSigma,1:2),intent(inout),optional :: sigma -- List of Xsection which is being updated in this routine First index: : channel Second index : 1=sum over runs( Xsection of each run), 2=sum over runs(( Xsection of each run)**2)
- logical,intent(inout),optional :: newInit
- integer,intent(inout),optional :: runNumber
- real, optional :: identifier -- Plotted in column 1 of the output file, if not present than we plot runNumber there

The optional input variable sigma returns the evaluated cross sections.

If newInit=.true., then the input sigma is initialized again (=0). If newInit=.false., then the input is used as a result of earlier runs and the results of the present run are just added.

If .true. then we assume that runNumber is the number to use as a normalization for the histogram results; We assume that the Xsection contain the results of "runNumber" runs.

NOTES

- If there are more than one particle of the same species produced, then the event is not considered for the "single particle cross section". e.g. a pi^+ pi^0 event does not contribute to any of those cross sections.
- If there are more than two particles of the same species produced, then the event is not considered for the "double particle cross section".
- If there are more than three particles of the same species produced, then the event is not considered for the "three particle cross section".
- If there are more than four particles of the same species produced, then the event is not considered for the "four particle cross section".

Note that e.g. 2pi/3pi means that two/three pions, no matter which charge, where found.

Note also that, e.g., "pi+ X" means that one piPlus was found. And if there is another particle of the same species, but with another charge, then we still count the event. E.g. for pi^+ pi^0 the event contributes both to "pi+ + X" and "pi0 + X"

OUTPUT

- filename 'sigma.dat'

## AnaEvent/sigma.dat [ Output files ]

[ Top ] [ AnaEvent ] [ Output files ]

NAME

file **sigma.dat**

PURPOSE

- First column : Naming for different runs or runNumber
- Columns 2-... : Single and multiparticle cross sections. The header of the file names the explicit channels.

NOTES

Cf. documenation in AnaEvent/event_sigma for the definition of "single" and "multi" particle cross sections. Esp. note the exclusive aspects.

## AnaEvent/event_dSigma_dOmega [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_dSigma_dOmega**(E, dTheta, dPhi, string, runNumber,
hists_theta, hists_phi, hists2D, initHists, makeOutputIn, sameFileNameIn)

PURPOSE

Output for single particle production dsigma/dTheta and dsigma/dTheta/dPhi with dTheta,dPhi in units of sr=radian for all particle-IDs includede in "particleIDs" (see above) and all charge states.

INPUTS

- type(tAnaEvent), dimension(:) :: E -- List of events
- real :: dTheta,dPhi -- Delta(Theta) and Delta(Phi) in [Radian]
- character(*) :: string -- used as prefix for all output file
- integer :: runNumber -- number of the run =1,2,3,4...

Optional:

- type(histogram) intent(inout), dimension(lBound(particleIDs,dim=1):uBound(particleIDs,dim=1),-2:2) :: hists_theta,hists_phi
- type(histogram2D), intent(inout), dimension(lBound(particleIDs,dim=1):uBound(particleIDs,dim=1),-2:2) :: hists2D
- logical, intent(in) :: initHists
- logical :: makeOutputIn -- Flag to switch off output (.false.=no output)
- logical, intent(in) :: sameFileNameIn -- .true. = always print to same filenames, .false. = different files for different runs

These optional input variables return the evaluated histograms. In hists_theta you will find the dsigma_dtheta histograms, in hists_phi you will find the dsigma_dphi histograms and in hists2D dsigma_dOmega. They must be given all together, to make an effect on the program.

If initHists=.true., then the input histograms "hists_theta","hists_phi" and "hists2D" are initialized again. If inithists=.false. then the input histograms are used as results of earlier runs and the results of the present run are just added. If .true. then we assume that runNumber is the number to use as a normalization for the histogram results = We assume that the histograms contain the results of "runNumber" runs.

The fourth column in the histogram defines the error of the 2nd column (only if optional input is given!)

NOTES

If there are more than one particle of the same species produced, then the event is not considered for this "single particle cross section".

OUTPUT

- filenames 'STRING_dsigma_dTheta...*.#.dat', 'STRING_dsigma_dTheta_Phi...*.#.dat' where #=1,2,3,... is given by runNumber and "STRING" by the input string.

## AnaEvent/event_dSigma_dE [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_dSigma_dE**(E, eMin, eMax, dE, string, runNumber, hists,
initHists, makeOutputIn, sameFileNameIn, histsMulti, hists1X, hists2X)

PURPOSE

Output for single particle production dsigma_dE for all particle IDs included in the array "particleIDs" (see above) and all charge states.

INPUTS

- type(tAnaEvent), dimension(:) :: E -- List of events
- real :: eMin, eMax, dE -- Minimal/Maximal Ekin, Delta(Ekin) in [GeV]
- character(*) :: string -- used as prefix for all output files
- integer :: runNumber -- number of the run

Optional:

- type(histogram), intent(inout), dimension(lBound(particleIDs,dim=1):uBound(particleIDs,dim=1),-2:2) :: hists
- logical, intent(in) :: initHists
- logical :: makeOutputIn -- Flag to switch off output (.false.=no output)
- logical, intent(in) :: sameFileNameIn -- .true. = always print to same filenames, .false. = different files for different runs

These optional input variables return the evaluated histograms. In hists you will find the dsigma_dE.

If initHists=.true., then the input histograms "hists" and "histsMulti" are initialized again. If inithists=.false., then the input histograms are used as results of earlier runs and the results of the present run are just added. If .true. then we assume that runNumber is the number to use as a normalization for the histogram results = We assume that the histograms contain the results of "runNumber" runs.

NOTES

If there are more than one particle of the same species produced, then the event is not considered for the "single particle cross section", but for the "multi particle cross section".

OUTPUT

- filenames 'STRING_dsigma_dEkin...*.#.dat' where #=1,2,3,... is given by runNumber and "STRING" by the input string.

## AnaEvent/event_dSigma_dEcostheta [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_dSigma_dEcostheta**(E, EcosthetaMin, EcosthetaMax,
dEcostheta, string, runNumber, hists, hists_c, hists_MULTI, hists_c_MULTI,
initHists, makeOutputIn, sameFileNameIn)

PURPOSE

Output for single particle production dsigma/d(E(1-cosTheta)) for all particle-IDs includede in "particleIDs" (see above) and all charge states. Theta is the angle of the outgoing particle with respect to the neutrino beam direction

The output is also given for dsigma/dcosTheta

INPUTS

- type(tAnaEvent), dimension(:) :: E -- List of events
- real ::
- character(*), intent(in) :: string -- used as prefix for all output file
- integer, intent(in) :: runNumber -- number of the run =1,2,3,4...

Optional:

- type(histogram), intent(inout), dimension(lBound(particleIDs,dim=1):uBound(particleIDs,dim=1),-2:2) :: hists
- logical, intent(in) :: initHists
- logical ::makeOutputIn -- Flag to switch off output (.false.=no output)
- logical , intent(in), optional ::sameFileNameIn -- .true. = always print to same filenames, .false. = different files for different runs

These optional input variables return the evaluated histograms.

If initHists=.true., then the input histograms are initialized again. If inithists=.false. then the input histograms are used as results of earlier runs and the results of the present run are just added. If .true. then we assume that runNumber is the number to use as a normalization for the histogram results = We assume that the histograms contain the results of "runNumber" runs.

The fourth column in the histogram defines the error of the 2nd column (only if optional input is given!)

NOTES

If there are more than one particle of the same species produced, then the event is not considered for this "single particle cross section".

OUTPUT

- filenames 'STRING_dsigma_dEcostheta...*.#.dat', where #=1,2,3,... is given by runNumber and "STRING" by the input string.

## AnaEvent/event_dSigma_dInvMass [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_dSigma_dInvMass**(E, string, runNumber,
dW_Npi, Wmin_Npi, Wmax_Npi, histsW_nucleon_pion,
dW_mupi, Wmin_mupi, Wmax_mupi, histsW_muon_pion,
dW_muN, Wmin_muN, Wmax_muN, histsW_muon_nucleon,
initHists, makeOutputIn, sameFileNameIn)

PURPOSE

Output dsigma_d(InvariantMass) for 1 pion production process for 3 invariant masses: W(pion-nucleon'), W(muon-nucleon'), W(pion-nucleon') can be compared with the corresponding ANL and BNL data

INPUTS

- type(tAnaEvent), dimension(:) :: E -- List of events
- real :: dInvMass -- Delta(InvMass) in [GeV]
- character(*) :: string -- used as prefix for all output files
- integer :: runNumber -- number of the run

Optional:

- type(histogram), intent(inout), dimension(lBound(particleIDs,dim=1):uBound(particleIDs,dim=1),-2:2) :: hists
- logical, intent(in) :: initHists
- logical :: makeOutputIn -- Flag to switch off output (.false.=no output)
- logical, intent(in) :: sameFileNameIn -- .true. = always print to same filenames, .false. = different files for different runs

These optional input variables return the evaluated histograms. In hists you will find the dsigma_dInvMass.

If initHists=.true., then the input histograms "hists" is initialized again. If inithists=.false., then the input histograms are used as results of earlier runs and the results of the present run are just added. If .true. then we assume that runNumber is the number to use as a normalization for the histogram results = We assume that the histograms contain the results of "runNumber" runs.

OUTPUT

- filenames 'STRING_dsigma_W...*.#.dat' where #=1,2,3,... is given by runNumber and "STRING" by the input string.

## AnaEvent/event_dSigma_dLeptonVariables [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_dSigma_dLeptonVariables**(E, eMin, eMax, dE, cost_min,
cost_max, delta_cost,string,
runNumber, hist_Enu, hist_E, hist_cos, hist_Q2, hist_Q2p, hist_dEdcost,hist_dpLdpT,
initHists, makeOutputIn, sameFileNameIn, specificEvent)

PURPOSE

Output dsi/dEkin, dsi/dcostheta for outgoing lepton and dsi/dQ^2

INPUTS

- type(tAnaEvent), dimension(:) :: E -- List of events
- real :: eMin, eMax, dE -- Minimal/Maximal Ekin, Delta(Ekin) in [GeV]
- real :: cost_min, cost_max, delta_cost -- limits for lepton scatt angle
- character(*) :: string -- used as prefix for all output files
- integer :: runNumber -- number of the run

Optional:

- type(histogram), intent(inout), dimension(lBound(particleIDs,dim=1):uBound(particleIDs,dim=1),-2:2) :: hist_Enu, hist_E, hist_cos, hist_Q2, hist_Q2p
- type(histogram2D),intent(inout) :: hist_dEdcost
- logical, intent(in) :: initHists
- logical :: makeOutputIn -- Flag to switch off output (.false.=no output)
- integer, intent(in) :: specificEvent -- specify final state (see subroutine SpecificEvent_Name)

These optional input variables return the evaluated histograms. In hist_E you will find the dsigma_dE(kinetic energy of the outgoing lepton). In hist_cos you will find the dsigma_dcostheta(angle of the outgoing lepton relative to incoming neutrino). In hist_Q2 you will find the dsigma_dQ^2. In hist_Q2p you will find the dsigma_dQ_p^2, Q_p^2 is calculated from kinetic energy of leading proton The Q_p^2 distribution can thus also be used to get the 0pi kinetic energy distribution of outgoing protons In hist_dEdcost you will find the double-differential d2sigma/(dE dcostheta) for outgoing lepton in a 0-pion event

If initHists=.true., then the input histograms "hists" and "histsMulti" are initialized again. If inithists=.false., then the input histograms are used as results of earlier runs and the results of the present run are just added. If .true. then we assume that runNumber is the number to use as a normalization for the histogram results = We assume that the histograms contain the results of "runNumber" runs.

OUTPUT

- filenames 'STRING_dsigma_dEkin_lepton.#.dat'
- filenames 'STRING_dsigma_dcos_lepton.#.dat'
- filenames 'STRING_dsigma_dQ2_lepton.#.dat'
- filenames 'STRING_dsigma_dQ2p_lepton.#.dat' where #=1,2,3,... is given by runNumber and "STRING" by the input string.

## AnaEvent/event_dump [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_dump**(run,L,filename,onlyFreeParticles,writeNeutrinoProd)

PURPOSE

Output a list of events.

INPUTS

- type(tAnaEvent),dimension(:), intent(in) :: L -- The eventlist
- character(100),intent(in) :: FileName -- The file name to write to
- integer,intent(in) :: run -- run number
- logical, optional :: onlyFreeParticles -- if .true. then only "free particles" are dumped. Here "free" means that p(0)^2-m_free^2 > 0. So far this is only implemented for nucleons.
- logical, optional :: writeNeutrinoProd -- switch to write production info

## AnaEvent/event_INIT [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_INIT**(L)

PURPOSE

Set the multiplicity counters back to zero.

INPUTS

- type(tAnaEvent) :: L -- The event

OUTPUT

(none)

## AnaEvent/event_CLEAR [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_CLEAR**(L)

PURPOSE

Reset the event: Delete the particleList and clear its memory. Set the multiplicity counters back to zero.

INPUTS

- type(tAnaEvent) :: L -- The event

OUTPUT

(none)

## AnaEvent/event_add [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_add**(L,V)

PURPOSE

Adds the particle (which V points at) to the event.

INPUTS

OUTPUT

(none)

## AnaEvent/event_getParticle [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

function **event_getParticle**(E,ID,Charge,n,P,antiparticle) RESULT (success)

PURPOSE

Search in the particle list of the event "E" the "n"-th particle in the list with %ID="ID". This particle "P" is then returned.

INPUTS

- type(tAnaEvent) :: E -- The event
- integer :: ID -- ID of particle which shall be returned
- integer :: Charge -- charge of particle which shall be returned
- integer :: n -- We return the n-th particle in the list with %ID=ID
- logical, optional :: antiparticle -- Whether the particle should be an antiparticle, if not specified, then we search for a particle.

OUTPUT

- type(particle) :: P -- n-th particle with wished ID
- logical :: success -- True if it was possible to find n-Particles with the wished ID, False otherwise

## AnaEvent/makeError_hist [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **makeError_hist**(a,b)

PURPOSE

Set b%y(:,3)=a%y(:,3)+( b%y(:,1)- a%y(:,1))**2

Mainly useful for error calculations. Usually performed to the histogram after a run with a being the histogram before the run and b being the histogram after the run.

INPUTS

- type(histogram), optional :: a -- Histogram of last run
- type(histogram), intent(inOut) :: b -- Histogram now

If a is missing then we assume that the histogram was 0 at the last run

OUTPUT

- type(histogram), intent(inOut) :: b -- Histogram now

## AnaEvent/setError_hist [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **setError_hist**(a,numRuns)

PURPOSE

Assume that a%yval(:,1) contains the following : sum over all runs(xSection of run). And assume that a%yval(:,3) contains the following : sum over all runs(xSection of run**2). The number of runs is given by numRuns.

Sets a%yval(:,3) to the statistical error of a%yval(:,1) times the number of runs (This is done because a%yval(:,1) is the mean value times the number of runs, too.

INPUTS

OUTPUT

## AnaEvent/event_DumpNumbers [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_DumpNumbers**

PURPOSE

...

INPUTS

- type(tAnaEvent) :: E ! The event

OUTPUT

- to stdout

## AnaEvent/event_pairPhotons [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_pairPhotons**(L,filename,dsigma_dm_threeGammas, &

& Hist_massPhotons)

PURPOSE

Outputs number of photons in one event and a histogram which gives the mass of all possible pairs of photons versus the mass of all possible triples. Additionally a cross section dsigma/dm_{3\gammas} is generated.

This is meant to reproduce an analysis procedure by the TAPS group.

INPUTS

- type(tAnaEvent),dimension(:) :: L -- The eventlist
- character(100) :: fileName_data,fileName_statistics -- The file name to write to

## AnaEvent/event_hadronicEnergy [ Subroutines ]

[ Top ] [ AnaEvent ] [ Subroutines ]

NAME

subroutine **event_hadronicEnergy**(L,filename,Ehad_versusNu, dSigdNu, dSigdEhad)

PURPOSE

Outputs histograms:

- Ehad_versusNu (2D) --- transferred energy (col#1); hadronic energy(col#2) xsec(col#3)
- dSigdNu --- dsi/dnu
- dSigdEhad --- dsi/dEhad
- Enurestored_versusEnu --- neutrino energy(col#1); restored neutrino energy (col#2) xsec(col#3)
- dSigdEnu --- dsi/dEnu
- dSigdEnurestored --- dsi/dEnurestored

Additionally info "prod_id Enu Enu_restored nu Ehad perweight" for all events is written into the file

This is meant to reproduce the calorimetric analysis by the MINOS group.

NOTES

The results are averaged over num_Energies (and, of course, also over num_Runs_sameEnergy)

INPUTS

- type(tAnaEvent),dimension(:) :: L -- The eventlist
- real :: numin, numax, nubin -- min, max and bin for the nu,Ehad-related histograms
- real :: Enumin, Enumax, Enubin -- min, max and bin for the Enu,Enurestored-related histograms

OUTPUT