diff --git a/CMakeLists.txt b/CMakeLists.txt
index bea638b..78fe008 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -1,83 +1,80 @@
# Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
################################################################################
# This file is part of NUISANCE.
#
# NUISANCE is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# NUISANCE is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with NUISANCE. If not, see .
################################################################################
cmake_minimum_required (VERSION 2.8 FATAL_ERROR)
#Use the compilers found in the path
find_program(CMAKE_C_COMPILER NAMES $ENV{CC} gcc PATHS ENV PATH NO_DEFAULT_PATH)
find_program(CMAKE_CXX_COMPILER NAMES $ENV{CXX} g++ PATHS ENV PATH NO_DEFAULT_PATH)
project(NUISANCE)
include(ExternalProject)
set (NUISANCE_VERSION_MAJOR 3)
set (NUISANCE_VERSION_MINOR 0)
set (NUISANCE_VERSION_REVISION 0)
set (NUISANCE_VERSION_STRING "v${NUISANCE_VERSION_MAJOR}r${NUISANCE_VERSION_MINOR}")
if(${NUISANCE_VERSION_REVISION} STRGREATER "0")
set (NUISANCE_VERSION_STRING "${NUISANCE_VERSION_STRING}p${NUISANCE_VERSION_REVISION}")
endif()
#Set this to TRUE to enable build debugging messages
set(BUILD_DEBUG_MSGS TRUE)
include(${CMAKE_SOURCE_DIR}/cmake/cmessage.cmake)
include(${CMAKE_SOURCE_DIR}/cmake/cacheVariables.cmake)
cmessage(STATUS "CMAKE_INSTALL_PREFIX: \"${CMAKE_INSTALL_PREFIX}\"")
cmessage(STATUS "CMAKE_BUILD_TYPE: \"${CMAKE_BUILD_TYPE}\"")
################################################################################
# Check Dependencies
################################################################################
################################## ROOT ######################################
include(${CMAKE_SOURCE_DIR}/cmake/ROOTSetup.cmake)
################################# InputHandler #################################
include(${CMAKE_SOURCE_DIR}/cmake/InputHandlerSetup.cmake)
################################# Pythia6/8 ###################################
include(${CMAKE_SOURCE_DIR}/cmake/pythia6Setup.cmake)
-################################## FHICLCPP ####################################
-include(${CMAKE_SOURCE_DIR}/cmake/fhiclcppSetup.cmake)
-
-#Need this to be at the front
-LIST(REVERSE EXTRA_CXX_FLAGS)
+#Want this before fhiclcpp which will add the install directory
LIST(APPEND EXTRA_CXX_FLAGS -I${CMAKE_SOURCE_DIR}/src)
-LIST(REVERSE EXTRA_CXX_FLAGS)
+################################## FHICLCPP ####################################
+include(${CMAKE_SOURCE_DIR}/cmake/fhiclcppSetup.cmake)
################################## COMPILER ####################################
include(${CMAKE_SOURCE_DIR}/cmake/c++CompilerSetup.cmake)
################################################################################
add_subdirectory(src)
configure_file(cmake/setup.sh.in
"${PROJECT_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/setup.sh" @ONLY)
install(FILES
"${PROJECT_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/setup.sh" DESTINATION
${CMAKE_INSTALL_PREFIX})
add_subdirectory(config)
add_subdirectory(data)
diff --git a/checksamples.sh b/checksamples.sh
index a0356af..edb23e5 100644
--- a/checksamples.sh
+++ b/checksamples.sh
@@ -1,29 +1,29 @@
for line in $(cat ./src/FCN/SampleList.cxx); #$(grep compare ../src/FCN/SampleList.cxx);
do
if [[ "$line" == *"compare"* ]];
then
parsed=${line/'!name.compare'/}
parsed=${parsed/'('/}
parsed=${parsed/')'/}
parsed=${parsed/'"'/}
parsed=${parsed/'("'/}
parsed=${parsed/'")'/}
echo $parsed
fi
if [[ "$line" == *"find("* ]];
then
parsed=${line/'!name.find'/}
parsed=${parsed/'('/}
parsed=${parsed/')'/}
parsed=${parsed/'"'/}
parsed=${parsed/'("'/}
parsed=${parsed/'")'/}
- echo "MultiSample: $line"
+ echo "MultIEventProcessor: $line"
fi
done;
diff --git a/cmake/NEUTSetup.cmake b/cmake/NEUTSetup.cmake
index 619d09f..2f7c60b 100644
--- a/cmake/NEUTSetup.cmake
+++ b/cmake/NEUTSetup.cmake
@@ -1,114 +1,130 @@
# Copyright 2016 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
################################################################################
# This file is part of NUISANCE.
#
# NUISANCE is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# NUISANCE is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with NUISANCE. If not, see .
################################################################################
if(NEUT_ROOT STREQUAL "")
cmessage(FATAL_ERROR "Variable NEUT_ROOT is not defined. Please export environment variable NEUT_ROOT or configure with -DNEUT_ROOT=/path/to/NEUT. This must be set to point to a prebuilt NEUT instance.")
endif()
if(CERN STREQUAL "")
cmessage(FATAL_ERROR "Variable CERN is not defined. Please export environment variable CERN or configure with -DCERN=/path/to/CERNLIB. This must be set to point to a prebuilt CERNLIB instance.")
endif()
if(CERN_LEVEL STREQUAL "")
cmessage(FATAL_ERROR "Variable CERN_LEVEL is not defined. Please export environment variable CERN_LEVEL or configure with -DCERN_LEVEL=XXXX (likely to be 2005).")
endif()
-if(NOT IS_NEUT_54)
+if(${NEUT_VERSION} VERSION_LESS 5.4.0)
set(NEUT_LIB_DIR ${NEUT_ROOT}/lib/Linux_pc)
else()
set(NEUT_LIB_DIR ${NEUT_ROOT}/lib)
endif()
+
set(NEUT_CLASS ${NEUT_ROOT}/src/neutclass)
-LIST(APPEND EXTRA_CXX_FLAGS -DNEUT_ENABLED )
+LIST(APPEND EXTRA_CXX_FLAGS -DNEUT_ENABLED -DNEUT_VERSION=${NEUT_VERSION})
+
+if(${NEUT_VERSION} VERSION_GREATER 5.4.0)
+ LIST(APPEND EXTRA_CXX_FLAGS -DNEUT_COMMON_QEAV)
+endif()
LIST(APPEND EXTRA_CXX_FLAGS
-I${NEUT_ROOT}/include
-I${NEUT_ROOT}/src/neutclass
-I${NEUT_ROOT}/src/reweight)
LIST(APPEND EXTRA_LINK_DIRS
${NEUT_LIB_DIR}
${CERN}/${CERN_LEVEL}/lib
${NEUT_ROOT}/src/reweight)
-if(IS_NEUT_54)
+LIST(APPEND EXTRA_LIBS NReWeight)
+
+if(${NEUT_VERSION} VERSION_EQUAL 5.4.2)
+ LIST(APPEND EXTRA_LIBS
+ neutcore_5.4.2
+ nuccorspl_5.4.2 #typo in NEUT, may hopefully disappear
+ nuceff_5.4.2
+ partnuck_5.4.2
+ skmcsvc_5.4.2
+ tauola_5.4.2
+ HT2p2h_5.4.0
+ N1p1h_5.4.0)
+ LIST(APPEND EXTRA_CXX_FLAGS -DNEUT_COMMON_QEAV)
+elseif(${NEUT_VERSION} VERSION_EQUAL 5.4.0)
LIST(APPEND EXTRA_LIBS
NReWeight
neutcore_5.4.0
nuccorspl_5.4.0 #typo in NEUT, may hopefully disappear
nuceff_5.4.0
partnuck_5.4.0
skmcsvc_5.4.0
tauola_5.4.0
HT2p2h_5.4.0
- N1p1h_5.4.0
- jetset74
- pdflib804
- mathlib
- packlib
- pawlib)
+ N1p1h_5.4.0)
else()
LIST(APPEND EXTRA_LIBS
NReWeight
neutcore
nuccorrspl
nuceff
partnuck
skmcsvc
- tauola
- jetset74
- pdflib804
- mathlib
- packlib
- pawlib)
+ tauola)
endif()
+LIST(APPEND EXTRA_LIBS
+ jetset74
+ pdflib804
+ mathlib
+ packlib
+ pawlib
+ gfortran)
+
set(NEUT_ROOT_LIBS)
LIST(APPEND NEUT_ROOT_LIBS
${NEUT_CLASS}/neutctrl.so
${NEUT_CLASS}/neutfsivert.so)
# Check for new versions of NEUT with NUCLEON FSI
if(EXISTS "${NEUT_CLASS}/neutnucfsistep.so")
set(NEUT_NUCFSI 1)
LIST(APPEND EXTRA_CXX_FLAGS -DNEUT_NUCFSI_ENABLED)
LIST(APPEND NEUT_ROOT_LIBS
${NEUT_CLASS}/neutnucfsistep.so
${NEUT_CLASS}/neutnucfsivert.so
)
endif()
-if(NOT IS_NEUT_54)
+if(${NEUT_VERSION} VERSION_LESS 5.4.0)
LIST(APPEND NEUT_ROOT_LIBS
${NEUT_CLASS}/neutrootTreeSingleton.so)
endif()
LIST(APPEND NEUT_ROOT_LIBS
${NEUT_CLASS}/neutvtx.so
${NEUT_CLASS}/neutfsipart.so
${NEUT_CLASS}/neutpart.so
${NEUT_CLASS}/neutvect.so
)
foreach(OBJ ${NEUT_ROOT_LIBS})
LIST(APPEND EXTRA_SHAREDOBJS ${OBJ})
endforeach()
diff --git a/cmake/c++CompilerSetup.cmake b/cmake/c++CompilerSetup.cmake
index 0443100..a6dd6cc 100644
--- a/cmake/c++CompilerSetup.cmake
+++ b/cmake/c++CompilerSetup.cmake
@@ -1,83 +1,83 @@
# Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
################################################################################
# This file is part of NUISANCE.
#
# NUISANCE is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# NUISANCE is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with NUISANCE. If not, see .
################################################################################
set(CXX_WARNINGS -Wall )
-LIST(APPEND EXTRA_CXX_FLAGS -std=c++1y -fPIC "-D__FILENAME__='\"$(subst ${CMAKE_SOURCE_DIR}/,,$(abspath $<))\"'")
+LIST(APPEND EXTRA_CXX_FLAGS -Wall -Wextra -Werror -Wno-delete-non-virtual-dtor -Wno-unused "-D__FILENAME__='\"$(subst ${CMAKE_SOURCE_DIR}/,,$(abspath $<))\"'")
cmessage(DEBUG "EXTRA_CXX_FLAGS: ${EXTRA_CXX_FLAGS}")
string(REPLACE ";" " " STR_EXTRA_CXX_FLAGS "${EXTRA_CXX_FLAGS}")
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${STR_EXTRA_CXX_FLAGS} ${CXX_WARNINGS}")
set(CMAKE_CXX_FLAGS_DEBUG "${CMAKE_CXX_FLAGS_DEBUG} -O0")
set(CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE} -O3")
SET(STR_EXTRA_LINK_DIRS)
if(NOT EXTRA_LINK_DIRS STREQUAL "")
string(REPLACE ";" " -L" STR_EXTRA_LINK_DIRS "-L${EXTRA_LINK_DIRS}")
endif()
LIST(APPEND EXTRA_LIBS dl)
SET(STR_EXTRA_LIBS)
if(NOT EXTRA_LIBS STREQUAL "")
SET(STR_EXTRA_LIBS_NO_SCRUB_LINKOPTS)
string(REPLACE ";" " -l" STR_EXTRA_LIBS_NO_SCRUB_LINKOPTS "-l${EXTRA_LIBS}")
string(REPLACE "-l-" "-" STR_EXTRA_LIBS ${STR_EXTRA_LIBS_NO_SCRUB_LINKOPTS})
endif()
SET(STR_EXTRA_SHAREDOBJS)
if(NOT EXTRA_SHAREDOBJS STREQUAL "")
string(REPLACE ";" " " STR_EXTRA_SHAREDOBJS "${EXTRA_SHAREDOBJS}")
endif()
LIST(APPEND EXTRA_LINK_FLAGS -Wl,--no-as-needed)
SET(STR_EXTRA_LINK_FLAGS)
if(NOT EXTRA_LINK_FLAGS STREQUAL "")
string(REPLACE ";" " " STR_EXTRA_LINK_FLAGS "${EXTRA_LINK_FLAGS}")
endif()
cmessage(DEBUG "EXTRA_LINK_DIRS: ${STR_EXTRA_LINK_DIRS}")
cmessage(DEBUG "EXTRA_LIBS: ${STR_EXTRA_LIBS}")
cmessage(DEBUG "EXTRA_SHAREDOBJS: ${STR_EXTRA_SHAREDOBJS}")
cmessage(DEBUG "EXTRA_LINK_FLAGS: ${STR_EXTRA_LINK_FLAGS}")
if(NOT STR_EXTRA_LINK_DIRS STREQUAL "" AND NOT STR_EXTRA_LIBS STREQUAL "")
SET(CMAKE_DEPENDLIB_FLAGS "${STR_EXTRA_LINK_DIRS} ${STR_EXTRA_LIBS}")
endif()
if(NOT STR_EXTRA_SHAREDOBJS STREQUAL "")
SET(CMAKE_DEPENDLIB_FLAGS "${CMAKE_DEPENDLIB_FLAGS} ${STR_EXTRA_SHAREDOBJS}")
endif()
if(NOT EXTRA_LINK_FLAGS STREQUAL "")
if(NOT CMAKE_LINK_FLAGS STREQUAL "")
SET(CMAKE_LINK_FLAGS "${CMAKE_LINK_FLAGS} ${STR_EXTRA_LINK_FLAGS}")
else()
SET(CMAKE_LINK_FLAGS "${STR_EXTRA_LINK_FLAGS}")
endif()
endif()
if (VERBOSE)
cmessage (STATUS "C++ Compiler : ${CXX_COMPILER_NAME}")
cmessage (STATUS " flags : ${CMAKE_CXX_FLAGS}")
cmessage (STATUS " Release flags : ${CMAKE_CXX_FLAGS_RELEASE}")
cmessage (STATUS " Debug flags : ${CMAKE_CXX_FLAGS_DEBUG}")
cmessage (STATUS " Link Flags : ${CMAKE_LINK_FLAGS}")
cmessage (STATUS " Lib Flags : ${CMAKE_DEPENDLIB_FLAGS}")
endif()
diff --git a/cmake/cacheVariables.cmake b/cmake/cacheVariables.cmake
index 11f3ae4..0c45699 100644
--- a/cmake/cacheVariables.cmake
+++ b/cmake/cacheVariables.cmake
@@ -1,99 +1,99 @@
# Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
################################################################################
# This file is part of NUISANCE.
#
# NUISANCE is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# NUISANCE is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with NUISANCE. If not, see .
################################################################################
function(CheckAndSetDefaultEnv VARNAME DEFAULT CACHETYPE DOCSTRING ENVNAME)
#cmessage(DEBUG "Trying to assign variable ${VARNAME} into the cache.")
if(NOT DEFINED ${VARNAME})
if(DEFINED ENV{${ENVNAME}} AND NOT $ENV{${ENVNAME}} STREQUAL "")
set(${VARNAME} $ENV{${ENVNAME}} CACHE ${CACHETYPE} ${DOCSTRING})
cmessage(DEBUG " Read ${VARNAME} from ENVVAR ${ENVNAME} as $ENV{${ENVNAME}}.")
else()
set(${VARNAME} ${DEFAULT} CACHE ${CACHETYPE} ${DOCSTRING})
endif()
else()
set(${VARNAME} ${${VARNAME}} CACHE ${CACHETYPE} ${DOCSTRING})
unset(${VARNAME})
endif()
cmessage(CACHE "--Set cache variable: \"${VARNAME}\" to \"${${VARNAME}}\", in cache ${CACHETYPE}.")
endfunction()
function(CheckAndSetDefaultCache VARNAME DEFAULT CACHETYPE DOCSTRING)
# cmessage(DEBUG "Trying to assign variable ${VARNAME} into the cache.")
if(NOT DEFINED ${VARNAME})
set(${VARNAME} ${DEFAULT} CACHE ${CACHETYPE} ${DOCSTRING})
else()
set(${VARNAME} ${${VARNAME}} CACHE ${CACHETYPE} ${DOCSTRING})
unset(${VARNAME})
endif()
cmessage(CACHE "--Set cache variable: \"${VARNAME}\" to \"${${VARNAME}}\", in cache ${CACHETYPE}.")
endfunction()
function(CheckAndSetDefault VARNAME DEFAULT)
# cmessage(DEBUG "Trying to assign variable ${VARNAME}.")
if(NOT DEFINED ${VARNAME})
set(${VARNAME} ${DEFAULT} PARENT_SCOPE)
set(${VARNAME} ${DEFAULT})
endif()
cmessage(CACHE "--Set variable: \"${VARNAME}\" to \"${${VARNAME}}\".")
endfunction()
CheckAndSetDefaultCache(VERBOSE TRUE BOOL "Whether to configure loudly.")
set (CMAKE_SKIP_BUILD_RPATH TRUE)
#Changes default install path to be a subdirectory of the build dir.
#Can set build dir at configure time with -DCMAKE_INSTALL_PREFIX=/install/path
if(CMAKE_INSTALL_PREFIX STREQUAL "" OR CMAKE_INSTALL_PREFIX STREQUAL
"/usr/local")
set(CMAKE_INSTALL_PREFIX "${CMAKE_BINARY_DIR}/${CMAKE_SYSTEM_NAME}")
elseif(NOT DEFINED CMAKE_INSTALL_PREFIX)
set(CMAKE_INSTALL_PREFIX "${CMAKE_BINARY_DIR}/${CMAKE_SYSTEM_NAME}")
endif()
if(CMAKE_BUILD_TYPE STREQUAL "")
set(CMAKE_BUILD_TYPE DEBUG)
elseif(NOT DEFINED CMAKE_BUILD_TYPE)
set(CMAKE_BUILD_TYPE DEBUG)
endif()
# Misc
CheckAndSetDefaultCache(EXTRA_SETUP_SCRIPT "" PATH "The path to an extra script to inject into the NUISANCE setup script. <>")
CheckAndSetDefaultCache(USE_MINIMIZER TRUE INTERNAL "Whether we are using the ROOT minimization libraries. ")
CheckAndSetDefaultCache(USE_ROOT6 FALSE INTERNAL "Whether we are using the ROOT 6. ")
# NuWro
CheckAndSetDefaultCache(USE_NuWro FALSE BOOL "Whether to enable NuWro support. ")
CheckAndSetDefaultEnv(NUWRO "" PATH "Path to NuWro source tree root directory. Overrides environment variable \$NUWRO <>" NUWRO)
CheckAndSetDefaultEnv(NUWRO_INC "" PATH "Path to NuWro installed includes directory, needs to contain \"params_all.h\". Overrides environment variable \$NUWRO_INC <>" NUWRO_INC)
# NEUT
CheckAndSetDefaultCache(USE_NEUT FALSE BOOL "Whether to enable NEUT (reweight) support. Requires external libraries. ")
-CheckAndSetDefaultCache(IS_NEUT_54 FALSE BOOL "Whether to enabled NEUT is version 5.4 or greater. ")
+CheckAndSetDefaultEnv(NEUT_VERSION FALSE STRING "NEUT version string, e.g. 5.4.0. <5.4.0>" NEUT_VERSION)
CheckAndSetDefaultEnv(NEUT_ROOT "" PATH "Path to NEUT source tree root directory. Overrides environment variable \$NEUT_ROOT <>" NEUT_ROOT)
CheckAndSetDefaultEnv(CERN "" PATH "Path to CERNLIB source tree root directory that NEUT was built against. Overrides environment variable \$CERN <>" CERN)
CheckAndSetDefaultEnv(CERN_LEVEL "" STRING "CERNLIB Library version. Overrides environment variable \$CERN_LEVEL <>" CERN_LEVEL)
# Pythia
CheckAndSetDefaultEnv(PYTHIA6 "" PATH "Path to directory containing libPythia6.so. Overrides environment variable \$PYTHIA6 <>" PYTHIA6)
CheckAndSetDefault(NEED_PYTHIA6 FALSE)
CheckAndSetDefault(NEED_ROOTEVEGEN FALSE)
CheckAndSetDefault(NEED_ROOTPYTHIA6 FALSE)
diff --git a/config/global/nuis.global.config.fcl.in b/config/global/nuis.global.config.fcl.in
index 5b964a4..2682d54 100644
--- a/config/global/nuis.global.config.fcl.in
+++ b/config/global/nuis.global.config.fcl.in
@@ -1,18 +1,20 @@
plugins: {
search_paths: {
#Default plugin install directory
installed_plugins: [ @CMAKE_INSTALL_PREFIX@/plugins ]
}
}
data_dir: @CMAKE_INSTALL_PREFIX@/data
persistency: {
- default_output_file: nuis.out.root
+ default: {
+ output_file: nuis.out.root
+ }
}
global: {
sample: {
verbosity_default: Reticent
}
}
diff --git a/data/CMakeLists.txt b/data/CMakeLists.txt
new file mode 100644
index 0000000..aafcffb
--- /dev/null
+++ b/data/CMakeLists.txt
@@ -0,0 +1,2 @@
+SET(DATA_INSTALL_DIR data)
+add_subdirectory(nuA)
diff --git a/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRD16_3103.root b/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRD16_3103.root
new file mode 100644
index 0000000..e9af236
Binary files /dev/null and b/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRD16_3103.root differ
diff --git a/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRD26_537.root b/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRD26_537.root
new file mode 100644
index 0000000..965ec79
Binary files /dev/null and b/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRD26_537.root differ
diff --git a/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRL31_844.root b/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRL31_844.root
new file mode 100644
index 0000000..11e20a0
Binary files /dev/null and b/data/nuA/BubbleChamber/ANL/CCQE/ANL_CCQE_Data_PRL31_844.root differ
diff --git a/data/nuA/BubbleChamber/ANL/CCQE/CMakeLists.txt b/data/nuA/BubbleChamber/ANL/CCQE/CMakeLists.txt
new file mode 100644
index 0000000..b71f174
--- /dev/null
+++ b/data/nuA/BubbleChamber/ANL/CCQE/CMakeLists.txt
@@ -0,0 +1,5 @@
+SET(DATA_INSTALL_DIR ${DATA_INSTALL_DIR}/CCQE)
+
+FILE(GLOB DATA_FILES *.root)
+
+install(FILES ${DATA_FILES} DESTINATION ${DATA_INSTALL_DIR})
diff --git a/data/nuA/BubbleChamber/ANL/CMakeLists.txt b/data/nuA/BubbleChamber/ANL/CMakeLists.txt
new file mode 100644
index 0000000..0d377cb
--- /dev/null
+++ b/data/nuA/BubbleChamber/ANL/CMakeLists.txt
@@ -0,0 +1,3 @@
+SET(DATA_INSTALL_DIR ${DATA_INSTALL_DIR}/ANL)
+
+add_subdirectory(CCQE)
diff --git a/data/nuA/BubbleChamber/CMakeLists.txt b/data/nuA/BubbleChamber/CMakeLists.txt
new file mode 100644
index 0000000..e4bba08
--- /dev/null
+++ b/data/nuA/BubbleChamber/CMakeLists.txt
@@ -0,0 +1,3 @@
+SET(DATA_INSTALL_DIR ${DATA_INSTALL_DIR}/BubbleChamber)
+
+add_subdirectory(ANL)
diff --git a/data/nuA/CMakeLists.txt b/data/nuA/CMakeLists.txt
new file mode 100644
index 0000000..36a3bdc
--- /dev/null
+++ b/data/nuA/CMakeLists.txt
@@ -0,0 +1,4 @@
+SET(DATA_INSTALL_DIR ${DATA_INSTALL_DIR}/nuA)
+
+add_subdirectory(BubbleChamber)
+add_subdirectory(Nuclear)
diff --git a/data/nuA/Nuclear/CMakeLists.txt b/data/nuA/Nuclear/CMakeLists.txt
new file mode 100644
index 0000000..b7a435f
--- /dev/null
+++ b/data/nuA/Nuclear/CMakeLists.txt
@@ -0,0 +1,3 @@
+SET(DATA_INSTALL_DIR ${DATA_INSTALL_DIR}/Nuclear)
+
+add_subdirectory(T2K)
diff --git a/data/nuA/Nuclear/T2K/CC0Pi/CMakeLists.txt b/data/nuA/Nuclear/T2K/CC0Pi/CMakeLists.txt
new file mode 100644
index 0000000..18cc6d1
--- /dev/null
+++ b/data/nuA/Nuclear/T2K/CC0Pi/CMakeLists.txt
@@ -0,0 +1,5 @@
+SET(DATA_INSTALL_DIR ${DATA_INSTALL_DIR}/CC0Pi)
+
+FILE(GLOB DATA_FILES *.root)
+
+install(FILES ${DATA_FILES} DESTINATION ${DATA_INSTALL_DIR})
diff --git a/data/nuA/Nuclear/T2K/CC0Pi/H2O_xsec_2Dpmuthetamu_numubar.root b/data/nuA/Nuclear/T2K/CC0Pi/H2O_xsec_2Dpmuthetamu_numubar.root
new file mode 100644
index 0000000..cae42cc
Binary files /dev/null and b/data/nuA/Nuclear/T2K/CC0Pi/H2O_xsec_2Dpmuthetamu_numubar.root differ
diff --git a/data/nuA/Nuclear/T2K/CMakeLists.txt b/data/nuA/Nuclear/T2K/CMakeLists.txt
new file mode 100644
index 0000000..25096b5
--- /dev/null
+++ b/data/nuA/Nuclear/T2K/CMakeLists.txt
@@ -0,0 +1,3 @@
+SET(DATA_INSTALL_DIR ${DATA_INSTALL_DIR}/T2K)
+
+add_subdirectory(CC0Pi)
diff --git a/doc/design_choices.md b/doc/design_choices.md
index b688f6a..86ad1f7 100644
--- a/doc/design_choices.md
+++ b/doc/design_choices.md
@@ -1,30 +1,30 @@
# Design Choices
This section describes the guiding design choices that are attempted to be
implemented in the NUISANCE v3 rewrite. Any perceived failure to adhere to these design choices should be regarded as a bug and submitted to the developers.
## Modularity of components
The various moving parts that go into a NUISANCE analysis should be as independent as possible (and no more). For example, as in previous versions, 'input handlers' are implementations of an `InputHandler` ABC, they are responsible for converting a specific event format into the NUISANCE event format. This event format should be generator agnostic as far as possible, however, for the generator-based-event reweighting to perform efficiently, the original generator event must be associated to the NUISANCE event. Instead of having some per-generator event manager, or reweight engines being allowed specific access to known InputHandler subclasses, it was decided that a pointer to the generator event would be carried around by the NUISANCE event base-class. This design choice has been kept.
-In NUISANCE v3, the design of the 'sample', the basic implementation unit of a comparison, has been further modularized. The previous Measurement{Base,1D,2D} base classes have been surplanted by the ISample and IDataComparison interfaces. IDataComparison being a subclass of ISample. These specify the interface that any NUISANCE sample must implement. The implementation of any subclass is left up to the user, however, a helper base class that should be used for almost all simple data comparisons can be found in SimpleDataComparison. It attempts to provide a fully functional base class, for data comparisons that previously subclassed from Measurement1D or 2D, that needs minimal specializing for specific data samples---most importantly the signal definition, the kinematic projection, and optional-but-encouraged data set metadata. As the interface suggests, the responsibility of looping over events and determining event variations is now fully that of the ISample, rather than a finely tuned collaboration between the sample and some calling method. However, subclasses of SimpleDataComparison can still be fully implemented by providing only per-event methods as the base class handles the event looping.
+In NUISANCE v3, the design of the 'sample', the basic implementation unit of a comparison, has been further modularized. The previous Measurement{Base,1D,2D} base classes have been surplanted by the IEventProcessor and IDataComparison interfaces. IDataComparison being a subclass of IEventProcessor. These specify the interface that any NUISANCE sample must implement. The implementation of any subclass is left up to the user, however, a helper base class that should be used for almost all simple data comparisons can be found in SimpleDataComparison. It attempts to provide a fully functional base class, for data comparisons that previously subclassed from Measurement1D or 2D, that needs minimal specializing for specific data samples---most importantly the signal definition, the kinematic projection, and optional-but-encouraged data set metadata. As the interface suggests, the responsibility of looping over events and determining event variations is now fully that of the IEventProcessor, rather than a finely tuned collaboration between the sample and some calling method. However, subclasses of SimpleDataComparison can still be fully implemented by providing only per-event methods as the base class handles the event looping.
Related to this, it may be useful to extend the 'reweight' concept, which was assumed to be the only form of event response in previous versions to a more general event 'variation'. However, since the majority of variations will still be weight-based, and as weight-based variations can often be calculated and applied fully with just the `nuis::event::MinimalEvent` format, the InputHandler will expose convenience methods that allow the current event weight to be without direct interrogation of any IWeightProvider instances. However, because general event variations can be arbitrarily complex and can require signal definitions to be rechecked (requiring a `nuis::event::FullEvent`), no such convenience method exists and samples looking to make use of such variations must explicitly call the VariationManager.
## Extensibility without recompilation
In previous versions of NUISANCE, the addition of new studies by non-experts was somewhat involved. Even though well-documented in multiple tutorials the addition of new studies required multiple source code changes, for opaque reasons outside of the study implementation itself. It also required recompilation of the entire set of NUISANCE binaries, this makes sharing a NUISANCE release on a cluster problematic and limits the entry barrier for non-expert users. V3 of NUISANCE is built from the ground up to allow for use and extensibility by non-experts, this is primarily facilitated by the `Instantiate` template method that implements a plugin factory, somewhat inspired by the `art::make_tool` utility from the ART framework. While the architecture of c++ means that plugins compiled with different toolsets than the main NUISANCE binaries are unlikely to work (and thus checked and disallowed before attempting any c++ object instantiation), multiple users working on a cluster should now be able to have local analysis implementations that are dynamically instantiated at runtime and require no recompilation of a central NUISANCE install.
## Flexibility of configuration
Originally NUISANCE used an inflexible DSL for configuration, in v2, a more-flexible XML-based configuration file was used. The change increased the scope of specific extensions to the core configuration, but the configuration API was clunky. In v3, NUISANCE uses an implementation of the FHiCL-c++ language
bindings to provide both configuration parsing and validation, as well as programatic access to the runtime configuration. FHiCL is particularly adept at providing a naturally hierarchial configuration, whereby global defaults can be
specified but overridden at runtime by user-supplied configuration files without special handling of the in-memory global configuration document. The FHiCL language bindings also specify how to retrieve configuration elements with specified types, reducing the need for comprehensive runtime checks for missing or incomplete configuration.
## Exceptional circumstances
-Previously, when NUISANCE encountered exceptional circumstances an exception was thrown. This exception was usually unqualified and sometimes preceeded by an error message. This made debugging and reporting of runtime errors fiddly for non-experts. In v3, all exceptions should be sensibly named, derived from `nuis::nuis_except` and be supplied with a verbose message diagnosing the problem encountered.
+Previously, when NUISANCE encountered exceptional circumstances an exception was thrown. This exception was usually unqualified and sometimes preceded by an error message. This made debugging and reporting of runtime errors fiddly for non-experts. In v3, all exceptions should be sensibly named, derived from `nuis::nuis_except` and be supplied with a verbose message diagnosing the problem encountered.
## Documentation
-The inline documentation of the code in previous versions of NUISANCE has been acceptable. The barrier to entry for new users was significantly lowered by external tutorials and documents. For v3, the hope is that the interface documentation can be fully comprehensive, and a specific set of instructions for common and less-common tasks can be distributed with the code. In addition to this, applications and scripts that verbosely aid users in standard NUISANCE workflows should be provided, and all usage prompts should be kept informative and fully up-to-date.
+The inline documentation of the code in previous versions of NUISANCE has been acceptable. The barrier to entry for new users was significantly lowered by external tutorials and documents. For v3, the hope is that the interface documentation can be fully comprehensive, and a specific set of instructions for common and less-common tasks can be distributed with the code. In addition to this, applications and scripts that verbosely aid users in standard NUISANCE workflows should be provided, and all usage prompts should be kept informative and fully up-to-date. Plugins should be queryable for example configurations.
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index b98e8d5..b615e75 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -1,13 +1,15 @@
add_subdirectory(config)
add_subdirectory(event)
add_subdirectory(input)
add_subdirectory(plugins)
add_subdirectory(utility)
add_subdirectory(generator)
add_subdirectory(persistency)
+add_subdirectory(parameters)
+add_subdirectory(variation)
add_subdirectory(samples)
SET(INuADataComparisons_List ${INuADataComparisons_List} PARENT_SCOPE)
SET(INuADataComparisons_FHiCL ${INuADataComparisons_FHiCL} PARENT_SCOPE)
add_subdirectory(app)
diff --git a/src/app/CMakeLists.txt b/src/app/CMakeLists.txt
index 3f12000..f11b727 100644
--- a/src/app/CMakeLists.txt
+++ b/src/app/CMakeLists.txt
@@ -1,13 +1,13 @@
SET(APPS nuissamples nuiscomp nuisevsum)
foreach(a ${APPS})
add_executable(${a} ${a}.cxx)
- target_link_libraries(${a} nuis_event nuis_input nuis_utility nuis_config nuis_persistency)
+ target_link_libraries(${a} nuis_event nuis_input nuis_utility nuis_config nuis_persistency nuis_plugins nuis_params nuis_variation)
target_link_libraries(${a} ${CMAKE_DEPENDLIB_FLAGS})
if(NOT "${CMAKE_LINK_FLAGS}" STREQUAL "")
set_target_properties(${a} PROPERTIES LINK_FLAGS ${CMAKE_LINK_FLAGS})
endif()
install(TARGETS ${a} DESTINATION bin)
endforeach()
diff --git a/src/app/nuiscomp.cxx b/src/app/nuiscomp.cxx
index da163f3..c69028f 100644
--- a/src/app/nuiscomp.cxx
+++ b/src/app/nuiscomp.cxx
@@ -1,53 +1,57 @@
#include "config/GlobalConfiguration.hxx"
#include "input/IInputHandler.hxx"
#include "event/MinimalEvent.hxx"
#include "samples/IDataComparison.hxx"
#include "plugins/Instantiate.hxx"
#include "exception/exception.hxx"
+#include "persistency/ROOTOutput.hxx"
+
#include "fhiclcpp/make_ParameterSet.h"
#include
NEW_NUIS_EXCEPT(invalid_cli_arguments);
int main(int argc, char const *argv[]) {
nuis::config::EnsureConfigurationRead("nuis.global.config.fcl");
nuis::config::EnsureConfigurationRead("nuis.datacomparisons.fcl");
if (argc != 2) {
throw invalid_cli_arguments()
<< "[ERROR]: Expected to be passed a single FHiCL file name or "
"absolute or relative path. N.B. Files in the local directory must "
"be fully qualified like \"$ "
<< argv[0] << " ./myconf.fcl\".";
}
nuis::config::EnsureConfigurationRead(argv[1]);
size_t NMax = nuis::config::GetDocument().get(
"nmax", std::numeric_limits::max());
for (fhicl::ParameterSet const &samp_config :
nuis::config::GetDocument().get>(
"samples")) {
std::cout << "[INFO]: Reading sample: "
<< samp_config.get("name") << std::endl;
nuis::plugins::plugin_traits::unique_ptr_t sample =
nuis::plugins::Instantiate(
samp_config.get("name"));
sample->Initialize(samp_config);
sample->ProcessSample(NMax);
std::cout << "[INFO]:\t Sample GOF = " << sample->GetGOF() << " / "
<< sample->GetNDOGuess() << std::endl;
sample->Write();
}
+
+ nuis::persistency::CloseOpenTFiles();
}
diff --git a/src/app/nuisevsum.cxx b/src/app/nuisevsum.cxx
index c45f677..5461e22 100644
--- a/src/app/nuisevsum.cxx
+++ b/src/app/nuisevsum.cxx
@@ -1,79 +1,79 @@
#include "config/GlobalConfiguration.hxx"
#include "input/IInputHandler.hxx"
#include "event/MinimalEvent.hxx"
-#include "samples/ISample.hxx"
+#include "samples/IEventProcessor.hxx"
#include "plugins/Instantiate.hxx"
#include "exception/exception.hxx"
#include "fhiclcpp/make_ParameterSet.h"
#include "string_parsers/from_string.hxx"
#include
NEW_NUIS_EXCEPT(invalid_cli_arguments);
size_t NMax = std::numeric_limits::max();
std::string input_file;
std::string input_type;
void SayUsage(char const *argv[]) {
std::cout << "[USAGE]: " << argv[0]
<< "\n"
"\t-i : Input file passed to named "
"IInputHandler instance \n"
"\t-H : Name of IInputHandler subclass "
"capable of reading NUISANCE events from the argument of -i.\n"
"\t-n : Maximum number of events to "
"read. Will read entire input file by default.\n"
<< std::endl;
}
void handleOpts(int argc, char const *argv[]) {
int opt = 1;
while (opt < argc) {
if ((std::string(argv[opt]) == "-?") ||
(std::string(argv[opt]) == "--help")) {
SayUsage(argv);
exit(0);
} else if (std::string(argv[opt]) == "-i") {
input_file = argv[++opt];
} else if (std::string(argv[opt]) == "-H") {
input_type = argv[++opt];
} else if (std::string(argv[opt]) == "-n") {
NMax = fhicl::string_parsers::str2T(argv[++opt]);
} else {
std::cout << "[ERROR]: Unknown option: " << argv[opt] << std::endl;
SayUsage(argv);
exit(1);
}
opt++;
}
}
int main(int argc, char const *argv[]) {
nuis::config::EnsureConfigurationRead("nuis.global.config.fcl");
handleOpts(argc, argv);
if (!input_type.length() || !input_file.length()) {
SayUsage(argv);
throw invalid_cli_arguments()
<< "[ERROR]: Require both -i and -H cli options to be passed.";
}
fhicl::ParameterSet sample_config;
sample_config.put("input_type", input_type);
sample_config.put("file", input_file);
- nuis::plugins::plugin_traits::unique_ptr_t VerboseEventSummary =
- nuis::plugins::Instantiate("VerboseEventSummary");
+ nuis::plugins::plugin_traits::unique_ptr_t VerboseEventSummary =
+ nuis::plugins::Instantiate("VerboseEventSummary");
VerboseEventSummary->Initialize(sample_config);
VerboseEventSummary->ProcessSample(NMax);
}
diff --git a/src/app/nuissamples.cxx b/src/app/nuissamples.cxx
index 3fdf60d..d3beeb5 100644
--- a/src/app/nuissamples.cxx
+++ b/src/app/nuissamples.cxx
@@ -1,111 +1,115 @@
#include "config/GlobalConfiguration.hxx"
#include "input/IInputHandler.hxx"
#include "plugins/Instantiate.hxx"
#include "samples/IDataComparison.hxx"
#include "utility/StringUtility.hxx"
#include
#include
#include
#include
std::string search_term = ".*";
bool fStrictRegex = false;
std::string config_out_filename = "";
void SayUsage(char const *argv[]) {
std::cout << "[USAGE]: " << argv[0]
<< "\n"
"\t-s : Used to filter known IDataComparisons. "
"\n\t\tWild cards are added to either side of the search term.\n"
"\t-S : Used to filter known IDataComparisons. "
"\n\t\tThe exact passed term is used.\n"
"\t-o : Dump example sample configuration file\n"
"\t\tfor matching samples.\n"
<< std::endl;
}
void handleOpts(int argc, char const *argv[]) {
int opt = 1;
while (opt < argc) {
if ((std::string(argv[opt]) == "-?") ||
(std::string(argv[opt]) == "--help")) {
SayUsage(argv);
exit(0);
} else if (std::string(argv[opt]) == "-s") {
search_term = argv[++opt];
} else if (std::string(argv[opt]) == "-S") {
fStrictRegex = true;
search_term = argv[++opt];
} else if (std::string(argv[opt]) == "-o") {
config_out_filename = argv[++opt];
} else {
std::cout << "[ERROR]: Unknown option: " << argv[opt] << std::endl;
SayUsage(argv);
exit(1);
}
opt++;
}
}
NEW_NUIS_EXCEPT(invalid_output_file);
int main(int argc, char const *argv[]) {
nuis::config::EnsureConfigurationRead("nuis.global.config.fcl");
nuis::config::EnsureConfigurationRead("nuis.datacomparisons.fcl");
handleOpts(argc, argv);
std::regex rpattern(fStrictRegex ? search_term
: std::string(".*") + search_term + ".*");
std::vector example_sample_configs;
for (std::string const &comparison_set_key :
nuis::config::GetDocument()
.get("data_comparisons")
.get_names()) {
for (std::string const &sample_name :
nuis::config::GetDocument().get>(
std::string("data_comparisons.") + comparison_set_key)) {
if (!std::regex_match(sample_name, rpattern)) {
continue;
}
nuis::plugins::plugin_traits::unique_ptr_t sample =
nuis::plugins::Instantiate(sample_name);
std::cout << sample->Name() << std::endl;
std::cout << "\tJournal: " << sample->GetJournalReference() << std::endl;
std::cout << "\tTarget: " << sample->GetTargetMaterial() << std::endl;
std::cout << "\tFlux: " << sample->GetFluxDescription() << std::endl;
std::cout << "\tSignal: " << sample->GetSignalDescription() << std::endl;
std::cout << "\tDocs: \n"
<< nuis::utility::indent_apply_width(sample->GetDocumentation(),
10)
<< std::endl;
- std::cout << "Config: "
- << sample->GetExampleConfiguration().to_indented_string()
+ std::cout << "\tExample_Config: {\n"
+ << nuis::utility::indent_apply_width(
+ sample->GetExampleConfiguration().to_indented_string(),
+ 12)
+ << "\n\t}\n"
<< std::endl;
+ ;
example_sample_configs.push_back(sample->GetExampleConfiguration());
}
}
if (config_out_filename.length()) {
std::ofstream out_file(config_out_filename);
if (!out_file) {
throw invalid_output_file() << "[ERROR]: Failed to open output file: "
<< std::quoted(config_out_filename);
}
fhicl::ParameterSet example_config;
example_config.put("samples", example_sample_configs);
out_file << example_config.to_indented_string();
}
}
diff --git a/src/config/GlobalConfiguration.cxx b/src/config/GlobalConfiguration.cxx
index ec35a15..772b429 100644
--- a/src/config/GlobalConfiguration.cxx
+++ b/src/config/GlobalConfiguration.cxx
@@ -1,41 +1,41 @@
#include "config/GlobalConfiguration.hxx"
#include "fhiclcpp/ParameterSet.h"
#include "fhiclcpp/make_ParameterSet.h"
namespace nuis {
namespace config {
struct NamedDoc {
std::string name;
fhicl::ParameterSet document;
std::vector files_read;
};
std::vector Configurations;
NamedDoc &GetDocument_modifyable(std::string const &name) {
for (NamedDoc &doc : Configurations) {
if (doc.name == name) {
return doc;
}
}
- Configurations.push_back(NamedDoc{name, fhicl::ParameterSet()});
+ Configurations.push_back(NamedDoc{name, fhicl::ParameterSet(), {}});
return GetDocument_modifyable(name);
}
fhicl::ParameterSet const &GetDocument(std::string const &name) {
return GetDocument_modifyable(name).document;
}
bool EnsureConfigurationRead(std::string const &fhicl_file,
std::string const &name) {
NamedDoc &doc = GetDocument_modifyable(name);
if (std::count(doc.files_read.begin(), doc.files_read.end(), fhicl_file)) {
return false;
}
doc.document.splice(fhicl::make_ParameterSet(fhicl_file));
return true;
}
} // namespace config
} // namespace nuis
diff --git a/src/event/FullEvent.cxx b/src/event/FullEvent.cxx
index 8d77987..c058989 100644
--- a/src/event/FullEvent.cxx
+++ b/src/event/FullEvent.cxx
@@ -1,24 +1,54 @@
#include "event/FullEvent.hxx"
namespace nuis {
namespace event {
FullEvent::FullEvent() : MinimalEvent() {
for (size_t status_it = 0;
status_it < static_cast(Particle::Status_t::kNParticleStatus);
++status_it) {
ParticleStack.push_back({static_cast(status_it), {}});
}
}
FullEvent::FullEvent(FullEvent &&other)
: MinimalEvent(std::move(other)),
ParticleStack(std::move(other.ParticleStack)) {}
+FullEvent &FullEvent::operator=(FullEvent &&other) {
+ MinimalEvent::operator=(std::move(other));
+ ParticleStack = std::move(other.ParticleStack);
+ return *this;
+}
+
+FullEvent FullEvent::Clone() const {
+ FullEvent clone;
+ clone.MinimalEvent::operator=(MinimalEvent::Clone());
+ clone.ParticleStack = ParticleStack;
+ return clone;
+}
+
void FullEvent::ClearParticleStack() {
for (auto &status_stack : ParticleStack) {
status_stack.particles.clear();
}
}
-} // namespace core
+std::string FullEvent::to_string() const {
+ std::stringstream ss("");
+ ss << "Event: Interaction mode = " << mode << ", probe: { PDG: " << probe_pdg
+ << ", Energy: " << probe_E << " MeV }." << std::endl;
+ for (auto &status_stack : ParticleStack) {
+ ss << "\t[" << status_stack.status << "]" << std::endl;
+
+ for (Particle const &part : status_stack.particles) {
+ ss << "\t\t{ PDG: " << part.pdg << ", P3: [ " << part.P4[0] << ", "
+ << part.P4[1] << ", " << part.P4[2] << "], E: " << part.P4[3]
+ << ", M: " << part.P4.M() << " }" << std::endl;
+ }
+ }
+ ss << std::endl;
+ return ss.str();
+}
+
+} // namespace event
} // namespace nuis
diff --git a/src/event/FullEvent.hxx b/src/event/FullEvent.hxx
index 74e4662..6e18f1a 100644
--- a/src/event/FullEvent.hxx
+++ b/src/event/FullEvent.hxx
@@ -1,70 +1,58 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef EVENT_FULLEVENT_HXX_SEEN
#define EVENT_FULLEVENT_HXX_SEEN
#include "event/MinimalEvent.hxx"
#include "event/Particle.hxx"
#include "string_parsers/to_string.hxx"
#include
namespace nuis {
namespace event {
///\brief The full, internal event format.
class FullEvent : public MinimalEvent {
public:
struct StatusParticles {
Particle::Status_t status;
std::vector particles;
};
FullEvent();
FullEvent(FullEvent const &) = delete;
FullEvent(FullEvent &&);
+ FullEvent& operator=(FullEvent &&);
+
+ FullEvent Clone() const;
+
std::vector ParticleStack;
void ClearParticleStack();
- std::string to_string() const {
- std::stringstream ss("");
- ss << "Event: Interaction mode = " << mode
- << ", probe: { PDG: " << probe_pdg << ", Energy: " << probe_E
- << " MeV }." << std::endl;
- for (auto &status_stack : ParticleStack) {
- ss << "\t[" << status_stack.status << "]" << std::endl;
-
- for (Particle const &part : status_stack.particles) {
- ss << "\t\t{ PDG: " << part.pdg << ", P3: [ " << part.P4[0] << ", "
- << part.P4[1] << ", " << part.P4[2] << "], E: " << part.P4[3]
- << ", M: " << part.P4.M() << " }" << std::endl;
- }
- }
- ss << std::endl;
- return ss.str();
- }
+ std::string to_string() const;
};
} // namespace event
} // namespace nuis
#endif
diff --git a/src/event/MinimalEvent.cxx b/src/event/MinimalEvent.cxx
index 1e78ed6..875b99f 100644
--- a/src/event/MinimalEvent.cxx
+++ b/src/event/MinimalEvent.cxx
@@ -1,22 +1,49 @@
#include "event/MinimalEvent.hxx"
namespace nuis {
namespace event {
MinimalEvent::MinimalEvent()
: mode(Channel_t::kUndefined), probe_E(0), probe_pdg(0), XSecWeight(1),
RWWeight(1) {
#ifdef __NUWRO_ENABLED__
fNuWroEvent = nullptr;
#endif
}
MinimalEvent::MinimalEvent(MinimalEvent &&other)
: mode(other.mode), probe_E(other.probe_E), probe_pdg(other.probe_pdg),
XSecWeight(other.XSecWeight), RWWeight(other.RWWeight) {
#ifdef __NUWRO_ENABLED__
fNuWroEvent = other.fNuWroEvent;
other.fNuWroEvent = nullptr;
#endif
}
-} // namespace core
+
+MinimalEvent &MinimalEvent::operator=(MinimalEvent &&other) {
+ mode = other.mode;
+ probe_E = other.probe_E;
+ probe_pdg = other.probe_pdg;
+ XSecWeight = other.XSecWeight;
+ RWWeight = other.RWWeight;
+#ifdef __NUWRO_ENABLED__
+ fNuWroEvent = other.fNuWroEvent;
+ other.fNuWroEvent = nullptr;
+#endif
+ return *this;
+}
+
+MinimalEvent MinimalEvent::Clone() const {
+ MinimalEvent clone;
+ clone.mode = mode;
+ clone.probe_E = probe_E;
+ clone.probe_pdg = probe_pdg;
+ clone.XSecWeight = XSecWeight;
+ clone.RWWeight = RWWeight;
+#ifdef __NUWRO_ENABLED__
+ clone.fNuWroEvent = fNuWroEvent;
+#endif
+
+ return clone;
+}
+} // namespace event
} // namespace nuis
diff --git a/src/event/MinimalEvent.hxx b/src/event/MinimalEvent.hxx
index 328965d..83f36f1 100644
--- a/src/event/MinimalEvent.hxx
+++ b/src/event/MinimalEvent.hxx
@@ -1,73 +1,78 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef EVENT_MINIMALEVENT_HXX_SEEN
#define EVENT_MINIMALEVENT_HXX_SEEN
#ifdef NUWRO_ENABLED
#include "event1.h"
typedef ::event NuWroEvent;
#endif
#ifdef NEUT_ENABLED
#include "neutpart.h"
#include "neutvect.h"
#endif
#include "event/types.hxx"
namespace nuis {
namespace event {
///\brief The minimal event information needed to perform reweights.
///
/// Most often, event selections cannot be applied using this reduced format.
class MinimalEvent {
public:
MinimalEvent();
MinimalEvent(MinimalEvent const &) = delete;
MinimalEvent(MinimalEvent &&);
+ MinimalEvent &operator=(MinimalEvent &&);
+
+ /// Make a clone of this MinimalEvent
+ MinimalEvent Clone() const;
+
/// True interaction mode
Channel_t mode;
/// True probe energy
double probe_E;
/// True probe particle code
PDG_t probe_pdg;
PDG_t target_pdg;
/// Event-weight that can be used to scale to a cross-section prediction
double XSecWeight;
/// Event weight incurred from current reweight engine state.
double RWWeight;
#ifdef NUWRO_ENABLED
///\brief Pointer to Nuwro event
///
/// This will usually be tied to a TTree and so we are not responsible for
/// deleting it
NuWroEvent *fNuWroEvent;
#endif
#ifdef NEUT_ENABLED
NeutVect *fNeutVect;
#endif
};
} // namespace core
} // namespace nuis
#endif
diff --git a/src/event/Particle.hxx b/src/event/Particle.hxx
index 0bcad11..dc6f4c0 100644
--- a/src/event/Particle.hxx
+++ b/src/event/Particle.hxx
@@ -1,77 +1,77 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef EVENT_PARTICLE_HXX_SEEN
#define EVENT_PARTICLE_HXX_SEEN
#include "event/types.hxx"
#include "TLorentzVector.h"
namespace nuis {
namespace event {
class Particle {
public:
NEW_NUIS_EXCEPT(invalid_particle);
#define STATUS_LIST \
X(kNuclearLeaving, 0) \
X(kPrimaryInitialState, 1) \
X(kPrimaryFinalState, 2) \
X(kIntermediate, 3) \
X(kUnknown, 4) \
X(kBlocked, 5) \
X(kNParticleStatus, 6)
#define X(A, B) A = B,
enum class Status_t { STATUS_LIST };
#undef X
Particle();
Particle(Particle const &);
PDG_t pdg;
TLorentzVector P4;
bool operator!() const { return (pdg == std::numeric_limits::max()); }
double E() const { return P4.E(); }
double P() const { return P4.Vect().Mag(); }
TVector3 P3() const { return P4.Vect(); }
double M() const { return P4.M(); }
double Theta() const { return P4.Vect().Theta(); }
- double CosTheta() const { return P4.Vect().Theta(); }
+ double CosTheta() const { return P4.Vect().CosTheta(); }
};
} // namespace event
} // namespace nuis
#define X(A, B) \
case nuis::event::Particle::Status_t::A: { \
return os << #A; \
}
inline std::ostream &operator<<(std::ostream &os,
nuis::event::Particle::Status_t te) {
switch (te) { STATUS_LIST }
return os;
}
#undef X
#undef STATUS_LIST
#endif
diff --git a/src/generator/CMakeLists.txt b/src/generator/CMakeLists.txt
index 84e4aea..eda12fb 100644
--- a/src/generator/CMakeLists.txt
+++ b/src/generator/CMakeLists.txt
@@ -1,2 +1,3 @@
add_subdirectory(input)
add_subdirectory(utility)
+add_subdirectory(variation)
diff --git a/src/generator/input/NEUTInputHandler.cxx b/src/generator/input/NEUTInputHandler.cxx
index a4133b6..798299e 100644
--- a/src/generator/input/NEUTInputHandler.cxx
+++ b/src/generator/input/NEUTInputHandler.cxx
@@ -1,103 +1,119 @@
#include "generator/input/NEUTInputHandler.hxx"
#include "utility/InteractionChannelUtility.hxx"
#include "utility/PDGCodeUtility.hxx"
#include "utility/ROOTUtility.hxx"
#include "generator/utility/NEUTUtility.hxx"
#include "fhiclcpp/ParameterSet.h"
using namespace nuis::event;
using namespace nuis::utility;
using namespace nuis::neuttools;
-NEUTInputHandler::NEUTInputHandler() : fInputTree(nullptr) {}
+NEUTInputHandler::NEUTInputHandler() : fInputTreeFile(nullptr) {}
NEUTInputHandler::NEUTInputHandler(NEUTInputHandler &&other)
- : fInputTree(std::move(other.fInputTree)),
+ : fInputTreeFile(std::move(other.fInputTreeFile)),
fReaderEvent(std::move(other.fReaderEvent)) {}
void NEUTInputHandler::Initialize(fhicl::ParameterSet const &ps) {
- fInputTree = CheckGetTTree(ps.get("file"), "neuttree");
+ fInputTreeFile = CheckGetTTree(ps.get("file"), "neuttree");
fReaderEvent.fNeutVect = nullptr;
- fInputTree->tree->SetBranchAddress("vectorbranch", &fReaderEvent.fNeutVect);
+ fInputTreeFile->tree->SetBranchAddress("vectorbranch",
+ &fReaderEvent.fNeutVect);
fKeepIntermediates = ps.get("keep_intermediates", false);
}
+double NEUTInputHandler::GetXSecScaleFactor(std::pair const &) const {
+ // check input file for histograms
+ std::unique_ptr flux =
+ GetHistogramFromROOTFile(fInputTreeFile->file, "flux_numu");
+ std::unique_ptr evtrt =
+ GetHistogramFromROOTFile(fInputTreeFile->file, "evtrt_numu");
+
+ if (flux && evtrt) {
+ return evtrt->Integral() / (flux->Integral() * double(GetNEvents()));
+ }
+
+ // try and build them
+ throw;
+}
+
MinimalEvent const &NEUTInputHandler::GetMinimalEvent(ev_index_t idx) const {
if (idx >= GetNEvents()) {
throw IInputHandler::invalid_entry()
<< "[ERROR]: Attempted to get entry " << idx
<< " from an InputHandler with only " << GetNEvents();
}
- fInputTree->tree->GetEntry(idx);
+ fInputTreeFile->tree->GetEntry(idx);
fReaderEvent.mode = IntToChannel(fReaderEvent.fNeutVect->Mode);
size_t NPart = fReaderEvent.fNeutVect->Npart();
for (size_t part_it = 0; part_it < NPart; part_it++) {
NeutPart const &part = (*fReaderEvent.fNeutVect->PartInfo(part_it));
if ((part.fIsAlive == false) && (part.fStatus == -1) &&
IsNeutralLepton(part.fPID)) {
fReaderEvent.probe_E = part.fP.T();
fReaderEvent.probe_pdg = part.fPID;
break;
}
}
return fReaderEvent;
}
FullEvent const &NEUTInputHandler::GetFullEvent(ev_index_t idx) const {
(void)GetMinimalEvent(idx);
fReaderEvent.ClearParticleStack();
if (fReaderEvent.fNeutVect->Ibound) {
fReaderEvent.target_pdg = MakeNuclearPDG(fReaderEvent.fNeutVect->TargetA,
fReaderEvent.fNeutVect->TargetZ);
} else {
fReaderEvent.target_pdg = MakeNuclearPDG(1, 1);
}
size_t NPart = fReaderEvent.fNeutVect->Npart();
size_t NPrimary = fReaderEvent.fNeutVect->Nprimary();
for (size_t part_it = 0; part_it < NPart; part_it++) {
NeutPart const &part = (*fReaderEvent.fNeutVect->PartInfo(part_it));
Particle nuis_part;
nuis_part.pdg = part.fPID;
nuis_part.P4 = part.fP;
Particle::Status_t state = GetNeutParticleStatus(part, fReaderEvent.mode);
if (!fKeepIntermediates && (state == Particle::Status_t::kIntermediate)) {
continue;
}
size_t state_int = static_cast(state);
// Add status == 0 particles as pre-FSI particles until we find an
// intermediate state particle
if ((part_it < NPrimary) &&
(state != Particle::Status_t::kPrimaryInitialState)) {
fReaderEvent
.ParticleStack[static_cast(
Particle::Status_t::kPrimaryFinalState)]
.particles.push_back(nuis_part);
}
fReaderEvent.ParticleStack[state_int].particles.push_back(nuis_part);
}
return fReaderEvent;
}
size_t NEUTInputHandler::GetNEvents() const {
- return fInputTree->tree->GetEntries();
+ return fInputTreeFile->tree->GetEntries();
}
DECLARE_PLUGIN(IInputHandler, NEUTInputHandler);
diff --git a/src/generator/input/NEUTInputHandler.hxx b/src/generator/input/NEUTInputHandler.hxx
index 9f4d8f6..4d063ea 100644
--- a/src/generator/input/NEUTInputHandler.hxx
+++ b/src/generator/input/NEUTInputHandler.hxx
@@ -1,56 +1,58 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef GENERATOR_INPUT_NEUTINPUTHANDLER_HXX_SEEN
#define GENERATOR_INPUT_NEUTINPUTHANDLER_HXX_SEEN
#include "event/FullEvent.hxx"
#include "input/IInputHandler.hxx"
#include
namespace fhicl {
class ParameterSet;
}
namespace nuis {
namespace utility {
class TreeFile;
}
} // namespace nuis
class NEUTInputHandler : public IInputHandler {
- mutable std::unique_ptr fInputTree;
+ mutable std::unique_ptr fInputTreeFile;
mutable nuis::event::FullEvent fReaderEvent;
bool fKeepIntermediates;
public:
NEUTInputHandler();
NEUTInputHandler(NEUTInputHandler const &) = delete;
NEUTInputHandler(NEUTInputHandler &&);
void Initialize(fhicl::ParameterSet const &);
nuis::event::MinimalEvent const &GetMinimalEvent(ev_index_t idx) const;
nuis::event::FullEvent const &GetFullEvent(ev_index_t idx) const;
size_t GetNEvents() const;
+ double GetXSecScaleFactor(
+ std::pair const &EnuRange) const;
};
#endif
diff --git a/src/generator/input/NuWroInputHandler.cxx b/src/generator/input/NuWroInputHandler.cxx
index 2fd3d52..52f2a7e 100644
--- a/src/generator/input/NuWroInputHandler.cxx
+++ b/src/generator/input/NuWroInputHandler.cxx
@@ -1,115 +1,121 @@
#include "generator/input/NuWroInputHandler.hxx"
#include "generator/utility/NuWroUtility.hxx"
#include "utility/ROOTUtility.hxx"
#include "fhiclcpp/ParameterSet.h"
#include "particle.h"
-typedef ::particle NuWroParticle;
+using NuWroParticle = ::particle;
using namespace nuis::event;
using namespace nuis::utility;
using namespace nuis::nuwrotools;
-NuWroInputHandler::NuWroInputHandler() : fInputTree(nullptr) {}
+NuWroInputHandler::NuWroInputHandler() : fInputTree() {}
NuWroInputHandler::NuWroInputHandler(NuWroInputHandler &&other)
: fInputTree(std::move(other.fInputTree)),
fReaderEvent(std::move(other.fReaderEvent)) {}
void NuWroInputHandler::Initialize(fhicl::ParameterSet const &ps) {
fInputTree = CheckGetTTree(ps.get("file"), "treeout");
fReaderEvent.fNuWroEvent = nullptr;
- fInputTree->tree->SetBranchAddress("e", &fReaderEvent.fNuWroEvent);
+ fInputTree.tree->SetBranchAddress("e", &fReaderEvent.fNuWroEvent);
fKeepIntermediates = ps.get("keep_intermediates", false);
}
MinimalEvent const &NuWroInputHandler::GetMinimalEvent(ev_index_t idx) const {
if (idx >= GetNEvents()) {
throw IInputHandler::invalid_entry()
<< "[ERROR]: Attempted to get entry " << idx
<< " from an InputHandler with only " << GetNEvents();
}
- fInputTree->tree->GetEntry(idx);
+ fInputTree.tree->GetEntry(idx);
fReaderEvent.mode = NuWroEventChannel(*fReaderEvent.fNuWroEvent);
fReaderEvent.probe_E = fReaderEvent.fNuWroEvent->in[0].E();
fReaderEvent.probe_pdg = fReaderEvent.fNuWroEvent->in[0].pdg;
fReaderEvent.XSecWeight =
fReaderEvent.fNuWroEvent->weight / double(GetNEvents());
if (fWeightCache.size() <= idx) {
fWeightCache.push_back(fReaderEvent.XSecWeight);
}
return fReaderEvent;
}
FullEvent const &NuWroInputHandler::GetFullEvent(ev_index_t idx) const {
(void)GetMinimalEvent(idx);
fReaderEvent.ClearParticleStack();
for (size_t p_it = 0; p_it < fReaderEvent.fNuWroEvent->in.size(); ++p_it) {
NuWroParticle &part = fReaderEvent.fNuWroEvent->in[p_it];
Particle nuis_part;
nuis_part.pdg = part.pdg;
nuis_part.P4 = TLorentzVector(part[1], part[2], part[3], part[0]);
fReaderEvent
.ParticleStack[static_cast(
Particle::Status_t::kPrimaryInitialState)]
.particles.push_back(nuis_part);
}
for (size_t p_it = 0;
p_it < fKeepIntermediates && fReaderEvent.fNuWroEvent->out.size();
++p_it) {
NuWroParticle &part = fReaderEvent.fNuWroEvent->out[p_it];
Particle nuis_part;
nuis_part.pdg = part.pdg;
nuis_part.P4 = TLorentzVector(part[1], part[2], part[3], part[0]);
fReaderEvent
.ParticleStack[static_cast(
Particle::Status_t::kPrimaryFinalState)]
.particles.push_back(nuis_part);
}
for (size_t p_it = 0; (p_it < fReaderEvent.fNuWroEvent->post.size());
++p_it) {
NuWroParticle &part = fReaderEvent.fNuWroEvent->post[p_it];
Particle nuis_part;
nuis_part.pdg = part.pdg;
nuis_part.P4 = TLorentzVector(part[1], part[2], part[3], part[0]);
fReaderEvent
.ParticleStack[static_cast(Particle::Status_t::kNuclearLeaving)]
.particles.push_back(nuis_part);
}
return fReaderEvent;
}
double NuWroInputHandler::GetEventWeight(ev_index_t idx) const {
if (idx > fWeightCache.size()) {
throw weight_cache_miss()
<< "[ERROR]: Failed to get cached weight for event index: " << idx;
}
return fWeightCache[idx];
}
+double NuWroInputHandler::GetXSecScaleFactor(
+ std::pair const &EnuRange) const {
+ throw input_handler_feature_unimplemented()
+ << "[ERROR]: Flux cuts not yet implemented for NuWro input handler.";
+}
+
size_t NuWroInputHandler::GetNEvents() const {
- return fInputTree->tree->GetEntries();
+ return fInputTree.tree->GetEntries();
}
DECLARE_PLUGIN(IInputHandler, NuWroInputHandler);
diff --git a/src/generator/input/NuWroInputHandler.hxx b/src/generator/input/NuWroInputHandler.hxx
index bb050e3..a7fc22c 100644
--- a/src/generator/input/NuWroInputHandler.hxx
+++ b/src/generator/input/NuWroInputHandler.hxx
@@ -1,63 +1,59 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
-#ifndef GENERATOR_INPUT_NUWROINPUTHANDLER_HXX_SEEN
-#define GENERATOR_INPUT_NUWROINPUTHANDLER_HXX_SEEN
+#pragma once
#include "event/FullEvent.hxx"
#include "input/IInputHandler.hxx"
#include "exception/exception.hxx"
+#include "utility/ROOTUtility.hxx"
+
#include
namespace fhicl {
class ParameterSet;
}
-namespace nuis {
-namespace utility {
-class TreeFile;
-}
-} // namespace nuis
-
class NuWroInputHandler : public IInputHandler {
- mutable std::unique_ptr fInputTree;
+ mutable nuis::utility::TreeFile fInputTree;
mutable nuis::event::FullEvent fReaderEvent;
mutable std::vector fWeightCache;
bool fKeepIntermediates;
public:
NEW_NUIS_EXCEPT(weight_cache_miss);
NuWroInputHandler();
NuWroInputHandler(NuWroInputHandler const &) = delete;
NuWroInputHandler(NuWroInputHandler &&);
void Initialize(fhicl::ParameterSet const &);
nuis::event::MinimalEvent const &GetMinimalEvent(ev_index_t idx) const;
nuis::event::FullEvent const &GetFullEvent(ev_index_t idx) const;
double GetEventWeight(ev_index_t idx) const;
size_t GetNEvents() const;
-};
-#endif
+ double GetXSecScaleFactor(
+ std::pair const &EnuRange) const;
+};
diff --git a/src/generator/utility/CMakeLists.txt b/src/generator/utility/CMakeLists.txt
index 7b60771..437ae87 100644
--- a/src/generator/utility/CMakeLists.txt
+++ b/src/generator/utility/CMakeLists.txt
@@ -1,20 +1,20 @@
if(USE_NuWro)
LIST(APPEND GENERATOR_UTILS_IMPL NuWroUtility.cxx)
LIST(APPEND GENERATOR_UTILS_HDR NuWroUtility.hxx)
endif(USE_NuWro)
if(USE_NEUT)
LIST(APPEND GENERATOR_UTILS_IMPL NEUTUtility.cxx)
LIST(APPEND GENERATOR_UTILS_HDR NEUTUtility.hxx)
endif(USE_NEUT)
if(GENERATOR_UTILS_IMPL)
add_library(nuis_generator_utility SHARED ${GENERATOR_UTILS_IMPL})
target_link_libraries(nuis_generator_utility nuis_event)
- install(TARGETS nuis_generator_utility DESTINATION lib)
+ install(TARGETS nuis_generator_utility DESTINATION plugins)
endif()
if(GENERATOR_UTILS_HDR)
install(FILES ${GENERATOR_UTILS_HDR} DESTINATION include/generator/utility)
endif()
diff --git a/src/generator/variation/CMakeLists.txt b/src/generator/variation/CMakeLists.txt
index e69de29..6cd1d8d 100644
--- a/src/generator/variation/CMakeLists.txt
+++ b/src/generator/variation/CMakeLists.txt
@@ -0,0 +1,13 @@
+LIST(APPEND WEIGHT_ENGINES_LINK_LIBS nuis_event nuis_config nuis_params)
+
+if(USE_NEUT)
+ LIST(APPEND WEIGHT_ENGINES_IMPL NEUTWeightEngine.cxx)
+ LIST(APPEND WEIGHT_ENGINES_LINK_LIBS nuis_generator_utility)
+endif(USE_NEUT)
+
+if(WEIGHT_ENGINES_IMPL)
+ add_library(generator_weight_engines SHARED ${WEIGHT_ENGINES_IMPL})
+ target_link_libraries(generator_weight_engines ${WEIGHT_ENGINES_LINK_LIBS})
+
+ install(TARGETS generator_weight_engines DESTINATION plugins)
+endif()
diff --git a/src/generator/variation/FillNEUTCommons.hxx b/src/generator/variation/FillNEUTCommons.hxx
new file mode 100644
index 0000000..8240d31
--- /dev/null
+++ b/src/generator/variation/FillNEUTCommons.hxx
@@ -0,0 +1,159 @@
+#ifndef GENERATOR_VARIATION_FILLENEUTCOMMONS_HXX_SEEN
+#define GENERATOR_VARIATION_FILLENEUTCOMMONS_HXX_SEEN
+
+#include "NFortFns.h" // Contains all the NEUT common blocks
+
+namespace NEUTUtils {
+
+inline void FillNeutCommons(NeutVect *nvect) {
+ // WARNING: This has only been implemented for a neuttree and not GENIE
+ // This should be kept in sync with T2KNIWGUtils::GetNIWGEvent(TTree)
+
+ // NEUT version info. Can't get it to compile properly with this yet
+ // neutversion_.corev = nvect->COREVer;
+ // neutversion_.nucev = nvect->NUCEVer;
+ // neutversion_.nuccv = nvect->NUCCVer;
+
+ // Documentation: See nework.h
+ nework_.modene = nvect->Mode;
+ nework_.numne = nvect->Npart();
+
+#ifdef NEUT_COMMON_QEAV
+ nemdls_.mdlqeaf = nvect->QEAVForm;
+#else
+ nemdls_.mdlqeaf = nvect->QEVForm;
+#endif
+ nemdls_.mdlqe = nvect->QEModel;
+ nemdls_.mdlspi = nvect->SPIModel;
+ nemdls_.mdldis = nvect->DISModel;
+ nemdls_.mdlcoh = nvect->COHModel;
+ neutcoh_.necohepi = nvect->COHModel;
+
+ nemdls_.xmaqe = nvect->QEMA;
+ nemdls_.xmvqe = nvect->QEMV;
+ nemdls_.kapp = nvect->KAPPA;
+
+ // nemdls_.sccfv = SCCFVdef;
+ // nemdls_.sccfa = SCCFAdef;
+ // nemdls_.fpqe = FPQEdef;
+
+ nemdls_.xmaspi = nvect->SPIMA;
+ nemdls_.xmvspi = nvect->SPIMV;
+ nemdls_.xmares = nvect->RESMA;
+ nemdls_.xmvres = nvect->RESMV;
+
+ neut1pi_.xmanffres = nvect->SPIMA;
+ neut1pi_.xmvnffres = nvect->SPIMV;
+ neut1pi_.xmarsres = nvect->RESMA;
+ neut1pi_.xmvrsres = nvect->RESMV;
+ neut1pi_.neiff = nvect->SPIForm;
+ neut1pi_.nenrtype = nvect->SPINRType;
+ neut1pi_.rneca5i = nvect->SPICA5I;
+ neut1pi_.rnebgscl = nvect->SPIBGScale;
+
+ nemdls_.xmacoh = nvect->COHMA;
+ nemdls_.rad0nu = nvect->COHR0;
+ // nemdls_.fa1coh = nvect->COHA1err;
+ // nemdls_.fb1coh = nvect->COHb1err;
+
+ // neutdis_.nepdf = NEPDFdef;
+ // neutdis_.nebodek = NEBODEKdef;
+
+ neutcard_.nefrmflg = nvect->FrmFlg;
+ neutcard_.nepauflg = nvect->PauFlg;
+ neutcard_.nenefo16 = nvect->NefO16;
+ neutcard_.nemodflg = nvect->ModFlg;
+ // neutcard_.nenefmodl = 1;
+ // neutcard_.nenefmodh = 1;
+ // neutcard_.nenefkinh = 1;
+ // neutpiabs_.neabspiemit = 1;
+
+ nenupr_.iformlen = nvect->FormLen;
+
+ neutpiless_.ipilessdcy = nvect->IPilessDcy;
+ neutpiless_.rpilessdcy = nvect->RPilessDcy;
+
+ neutpiless_.ipilessdcy = nvect->IPilessDcy;
+ neutpiless_.rpilessdcy = nvect->RPilessDcy;
+
+ neffpr_.fefqe = nvect->NuceffFactorPIQE;
+ neffpr_.fefqeh = nvect->NuceffFactorPIQEH;
+ neffpr_.fefinel = nvect->NuceffFactorPIInel;
+ neffpr_.fefabs = nvect->NuceffFactorPIAbs;
+ neffpr_.fefcx = nvect->NuceffFactorPICX;
+ neffpr_.fefcxh = nvect->NuceffFactorPICXH;
+
+ neffpr_.fefcoh = nvect->NuceffFactorPICoh;
+ neffpr_.fefqehf = nvect->NuceffFactorPIQEHKin;
+ neffpr_.fefcxhf = nvect->NuceffFactorPICXKin;
+ neffpr_.fefcohf = nvect->NuceffFactorPIQELKin;
+
+ for (int i = 0; i < nework_.numne; i++) {
+ nework_.ipne[i] = nvect->PartInfo(i)->fPID;
+ nework_.pne[i][0] =
+ (float)nvect->PartInfo(i)->fP.X() / 1000; // VC(NE)WORK in M(G)eV
+ nework_.pne[i][1] =
+ (float)nvect->PartInfo(i)->fP.Y() / 1000; // VC(NE)WORK in M(G)eV
+ nework_.pne[i][2] =
+ (float)nvect->PartInfo(i)->fP.Z() / 1000; // VC(NE)WORK in M(G)eV
+ }
+ // fsihist.h
+
+ // neutroot fills a dummy object for events with no FSI to prevent memory leak
+ // when
+ // reading the TTree, so check for it here
+
+ if ((int)nvect->NfsiVert() ==
+ 1) { // An event with FSI must have at least two vertices
+ // if (nvect->NfsiPart()!=1 || nvect->Fsiprob!=-1)
+ // ERR(WRN) << "T2KNeutUtils::fill_neut_commons(TTree) NfsiPart!=1 or
+ // Fsiprob!=-1 when NfsiVert==1" << std::endl;
+
+ fsihist_.nvert = 0;
+ fsihist_.nvcvert = 0;
+ fsihist_.fsiprob = 1;
+ } else { // Real FSI event
+ fsihist_.nvert = (int)nvect->NfsiVert();
+ for (int ivert = 0; ivert < fsihist_.nvert; ivert++) {
+ fsihist_.iflgvert[ivert] = nvect->FsiVertInfo(ivert)->fVertID;
+ fsihist_.posvert[ivert][0] = (float)nvect->FsiVertInfo(ivert)->fPos.X();
+ fsihist_.posvert[ivert][1] = (float)nvect->FsiVertInfo(ivert)->fPos.Y();
+ fsihist_.posvert[ivert][2] = (float)nvect->FsiVertInfo(ivert)->fPos.Z();
+ }
+
+ fsihist_.nvcvert = nvect->NfsiPart();
+ for (int ip = 0; ip < fsihist_.nvcvert; ip++) {
+ fsihist_.abspvert[ip] = (float)nvect->FsiPartInfo(ip)->fMomLab;
+ fsihist_.abstpvert[ip] = (float)nvect->FsiPartInfo(ip)->fMomNuc;
+ fsihist_.ipvert[ip] = nvect->FsiPartInfo(ip)->fPID;
+ fsihist_.iverti[ip] = nvect->FsiPartInfo(ip)->fVertStart;
+ fsihist_.ivertf[ip] = nvect->FsiPartInfo(ip)->fVertEnd;
+ fsihist_.dirvert[ip][0] = (float)nvect->FsiPartInfo(ip)->fDir.X();
+ fsihist_.dirvert[ip][1] = (float)nvect->FsiPartInfo(ip)->fDir.Y();
+ fsihist_.dirvert[ip][2] = (float)nvect->FsiPartInfo(ip)->fDir.Z();
+ }
+ fsihist_.fsiprob = nvect->Fsiprob;
+ }
+
+ neutcrscom_.crsx = nvect->Crsx;
+ neutcrscom_.crsy = nvect->Crsy;
+ neutcrscom_.crsz = nvect->Crsz;
+ neutcrscom_.crsphi = nvect->Crsphi;
+ neutcrscom_.crsq2 = nvect->Crsq2;
+
+ neuttarget_.numbndn = nvect->TargetA - nvect->TargetZ;
+ neuttarget_.numbndp = nvect->TargetZ;
+ neuttarget_.numfrep = nvect->TargetH;
+ neuttarget_.numatom = nvect->TargetA;
+ posinnuc_.ibound = nvect->Ibound;
+
+ // put empty nucleon FSI history (since it is not saved in the NeutVect
+ // format)
+ // Comment out as NEUT does not have the necessary proton FSI information yet
+ // nucleonfsihist_.nfnvert = 0;
+ // nucleonfsihist_.nfnstep = 0;
+}
+
+} // namespace NEUTUtils
+
+#endif
diff --git a/src/generator/variation/NEUTReWeight.cxx b/src/generator/variation/NEUTReWeight.cxx
deleted file mode 100644
index 4e61a9c..0000000
--- a/src/generator/variation/NEUTReWeight.cxx
+++ /dev/null
@@ -1,3 +0,0 @@
-#include "generator/variation/NEUTReWeight.hxx"
-
-DECLARE_PLUGIN(IWeightProvider, NEUTReWeight);
diff --git a/src/generator/variation/NEUTWeightEngine.cxx b/src/generator/variation/NEUTWeightEngine.cxx
new file mode 100644
index 0000000..51c1977
--- /dev/null
+++ b/src/generator/variation/NEUTWeightEngine.cxx
@@ -0,0 +1,114 @@
+#include "generator/variation/NEUTWeightEngine.hxx"
+#include "generator/variation/FillNEUTCommons.hxx"
+
+#include "fhiclcpp/ParameterSet.h"
+
+#include "event/MinimalEvent.hxx"
+
+// NEUT Engine includes
+#include "NReWeight.h"
+#include "NReWeightCasc.h"
+#include "NReWeightNuXSecCCQE.h"
+#include "NReWeightNuXSecCCRES.h"
+#include "NReWeightNuXSecCOH.h"
+#include "NReWeightNuXSecDIS.h"
+#include "NReWeightNuXSecNC.h"
+#include "NReWeightNuXSecNCEL.h"
+#include "NReWeightNuXSecNCRES.h"
+#include "NReWeightNuXSecRES.h"
+#include "NReWeightNuclPiless.h"
+
+#include "NSystUncertainty.h"
+
+#include
+
+using namespace nuis;
+using namespace nuis::params;
+using namespace nuis::event;
+
+void NEUTWeightEngine::Initialize(fhicl::ParameterSet const &ps) {
+
+ fNeutRW = std::make_unique();
+ TDirectory *dir = gDirectory;
+ fNeutRW->AdoptWghtCalc("xsec_ccqe", new neut::rew::NReWeightNuXSecCCQE);
+ fNeutRW->AdoptWghtCalc("xsec_res", new neut::rew::NReWeightNuXSecRES);
+ fNeutRW->AdoptWghtCalc("xsec_ccres", new neut::rew::NReWeightNuXSecCCRES);
+ fNeutRW->AdoptWghtCalc("xsec_coh", new neut::rew::NReWeightNuXSecCOH);
+ fNeutRW->AdoptWghtCalc("xsec_dis", new neut::rew::NReWeightNuXSecDIS);
+ fNeutRW->AdoptWghtCalc("xsec_ncel", new neut::rew::NReWeightNuXSecNCEL);
+ fNeutRW->AdoptWghtCalc("xsec_nc", new neut::rew::NReWeightNuXSecNC);
+ fNeutRW->AdoptWghtCalc("xsec_ncres", new neut::rew::NReWeightNuXSecNCRES);
+ fNeutRW->AdoptWghtCalc("nucl_casc", new neut::rew::NReWeightCasc);
+ fNeutRW->AdoptWghtCalc("nucl_piless", new neut::rew::NReWeightNuclPiless);
+ dir->cd();
+
+ for (fhicl::ParameterSet param_ps :
+ ps.get>("parameters")) {
+
+ param_ps.put("type", GetName());
+
+ std::string const ¶m_name = param_ps.get("name");
+
+ NEUTSystParam nsp;
+
+ nsp.nsyst = neut::rew::NSyst::FromString(param_name);
+
+ if (nsp.nsyst == neut::rew::kNullSystematic) {
+ throw invalid_NEUT_syst_name() << "[ERROR]: NReWeight failed to parse "
+ << std::quoted(param_name) << " as a NEUT dial.";
+ }
+ fNeutRW->Systematics().Init(nsp.nsyst);
+
+ nsp.pid = ParameterManager::Get().EnsureParameterRegistered(param_ps);
+
+ fNEUTSysts.push_back(nsp);
+ }
+
+ Reconfigure();
+}
+void NEUTWeightEngine::Reconfigure() {
+
+ for (NEUTSystParam const &nsp : fNEUTSysts) {
+ double val = ParameterManager::Get().GetParameterValue(nsp.pid);
+ fNeutRW->Systematics().Set(nsp.nsyst, val);
+ }
+
+ fNeutRW->Reconfigure();
+}
+double NEUTWeightEngine::GetEventWeight(nuis::event::MinimalEvent const &ev) {
+ if (!ev.fNeutVect) {
+ return 1.0;
+ }
+ NEUTUtils::FillNeutCommons(ev.fNeutVect);
+ return fNeutRW->CalcWeight();
+}
+
+std::string NEUTWeightEngine::GetName() { return "NEUTWeightEngine"; }
+std::string NEUTWeightEngine::GetDocumentation() { return ""; }
+fhicl::ParameterSet NEUTWeightEngine::GetExampleConfiguration() {
+ fhicl::ParameterSet ps;
+
+ ps.put("weight_engine_name", GetName());
+ fhicl::ParameterSet dial_maqe;
+ dial_maqe.put("name", "MAQE");
+
+ dial_maqe.put("start", 0);
+ dial_maqe.put("min", -3);
+ dial_maqe.put("max", 3);
+ dial_maqe.put("step", 0.1);
+
+ fhicl::ParameterSet dial_mares;
+ dial_mares.put("name", "MARES");
+
+ dial_mares.put("start", 0);
+ dial_mares.put("min", -3);
+ dial_mares.put("max", 3);
+ dial_mares.put("step", 0.1);
+
+ ps.put>(
+ "parameters", std::vector{{dial_maqe, dial_mares}});
+
+ return ps;
+}
+
+DECLARE_PLUGIN(IWeightProvider, NEUTWeightEngine);
diff --git a/src/generator/variation/NEUTReWeight.hxx b/src/generator/variation/NEUTWeightEngine.hxx
similarity index 67%
rename from src/generator/variation/NEUTReWeight.hxx
rename to src/generator/variation/NEUTWeightEngine.hxx
index fd2c330..b7581c8 100644
--- a/src/generator/variation/NEUTReWeight.hxx
+++ b/src/generator/variation/NEUTWeightEngine.hxx
@@ -1,35 +1,61 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef GENERATOR_VARIATION_NEUTREWEIGHT_HXX_SEEN
#define GENERATOR_VARIATION_NEUTREWEIGHT_HXX_SEEN
#include "variation/IWeightProvider.hxx"
-class NEUTReWeight : public IWeightProvider {
+#include "exception/exception.hxx"
+
+#include "parameters/ParameterManager.hxx"
+
+#include "NSyst.h"
+
+#include
+#include
+
+namespace neut {
+namespace rew {
+class NReWeight;
+}
+} // namespace neut
+
+class NEUTWeightEngine : public IWeightProvider {
+
+ struct NEUTSystParam {
+ nuis::params::paramId_t pid;
+ neut::rew::NSyst_t nsyst;
+ };
+
+ std::vector fNEUTSysts;
+ std::unique_ptr fNeutRW;
+
public:
- double GetEventWeight(nuis::event::MinimalEvent const &);
+ NEW_NUIS_EXCEPT(invalid_NEUT_syst_name);
+
void Initialize(fhicl::ParameterSet const &);
- paramId_t GetParameterId(std::string const &);
- void SetParameterValue(paramId_t, double);
- bool ParametersVaried();
void Reconfigure();
+ double GetEventWeight(nuis::event::MinimalEvent const &);
+ std::string GetName();
+ std::string GetDocumentation();
+ fhicl::ParameterSet GetExampleConfiguration();
};
#endif
diff --git a/src/input/CMakeLists.txt b/src/input/CMakeLists.txt
index 5aad051..6036777 100644
--- a/src/input/CMakeLists.txt
+++ b/src/input/CMakeLists.txt
@@ -1,12 +1,12 @@
SET(input_implementation_files
InputManager.cxx)
SET(input_header_files
InputManager.hxx
IInputHandler.hxx)
add_library(nuis_input SHARED ${input_implementation_files})
-target_link_libraries(nuis_input)
+target_link_libraries(nuis_input nuis_plugins)
install(TARGETS nuis_input DESTINATION lib)
install(FILES ${input_header_files} DESTINATION include/input)
diff --git a/src/input/IInputHandler.hxx b/src/input/IInputHandler.hxx
index e558904..58ea6c0 100644
--- a/src/input/IInputHandler.hxx
+++ b/src/input/IInputHandler.hxx
@@ -1,113 +1,118 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
-#ifndef INPUT_IINPUTHANDLER_HXX_SEEN
-#define INPUT_IINPUTHANDLER_HXX_SEEN
+#pragma once
#include "plugins/traits.hxx"
#include "exception/exception.hxx"
#include
+#include
namespace fhicl {
class ParameterSet;
}
namespace nuis {
namespace event {
class MinimalEvent;
class FullEvent;
} // namespace event
} // namespace nuis
class IInputHandler {
public:
struct FullEvent_const_iterator
: public std::iterator<
std::input_iterator_tag, nuis::event::FullEvent const, size_t,
nuis::event::FullEvent const *, nuis::event::FullEvent const &> {
FullEvent_const_iterator(size_t _idx, IInputHandler const *_ih) {
idx = _idx;
ih = _ih;
}
FullEvent_const_iterator(FullEvent_const_iterator const &other) {
idx = other.idx;
ih = other.ih;
}
FullEvent_const_iterator operator=(FullEvent_const_iterator const &other) {
idx = other.idx;
ih = other.ih;
return (*this);
}
bool operator==(FullEvent_const_iterator const &other) {
return (idx == other.idx);
}
bool operator!=(FullEvent_const_iterator const &other) {
return !(*this == other);
}
nuis::event::FullEvent const &operator*() { return ih->GetFullEvent(idx); }
nuis::event::FullEvent const *operator->() {
return &ih->GetFullEvent(idx);
}
FullEvent_const_iterator operator++() {
idx++;
return *this;
}
FullEvent_const_iterator operator++(int) {
FullEvent_const_iterator tmp(*this);
idx++;
return tmp;
}
private:
size_t idx;
IInputHandler const *ih;
};
NEW_NUIS_EXCEPT(invalid_input_file);
NEW_NUIS_EXCEPT(invalid_entry);
+ NEW_NUIS_EXCEPT(input_handler_feature_unimplemented);
typedef size_t ev_index_t;
virtual void Initialize(fhicl::ParameterSet const &) = 0;
virtual nuis::event::MinimalEvent const &
GetMinimalEvent(ev_index_t idx) const = 0;
virtual nuis::event::FullEvent const &GetFullEvent(ev_index_t idx) const = 0;
virtual void RecalculateEventWeights(){};
- virtual double GetEventWeight(ev_index_t idx) const { return 1; };
+ virtual double GetEventWeight(ev_index_t) const { return 1; };
+
+ /// Allows samples that implement flux cuts to request an appropriate
+ virtual double GetXSecScaleFactor(
+ std::pair const &EnuRange = std::pair{
+ std::numeric_limits::max(),
+ std::numeric_limits::max()}) const = 0;
virtual size_t GetNEvents() const = 0;
FullEvent_const_iterator begin() const {
return FullEvent_const_iterator(0, this);
}
FullEvent_const_iterator end() const {
return FullEvent_const_iterator(GetNEvents(), this);
}
virtual ~IInputHandler() {}
};
DECLARE_PLUGIN_INTERFACE(IInputHandler);
-
-#endif
diff --git a/src/input/InputManager.cxx b/src/input/InputManager.cxx
index 8928d84..6e2dbca 100644
--- a/src/input/InputManager.cxx
+++ b/src/input/InputManager.cxx
@@ -1,65 +1,66 @@
#include "input/InputManager.hxx"
#include "plugins/Instantiate.hxx"
#include "fhiclcpp/ParameterSet.h"
#include
namespace nuis {
namespace input {
InputManager *InputManager::_global_inst = nullptr;
InputManager::NamedInputHandler::NamedInputHandler(
std::string const &file,
plugins::plugin_traits::unique_ptr_t &&IH) {
name = file;
handler = std::move(IH);
}
InputManager::InputManager() : Inputs() {}
InputManager &InputManager::Get() {
if (!_global_inst) {
_global_inst = new InputManager();
}
return *_global_inst;
}
InputManager::Input_id_t
InputManager::EnsureInputLoaded(fhicl::ParameterSet const &ps) {
std::string const &file_name = ps.get("file");
for (size_t i = 0; i < Inputs.size(); ++i) {
if (Inputs[i].name == file_name) {
return i;
}
}
+ Input_id_t iid = Inputs.size();
Inputs.emplace_back(file_name,
nuis::plugins::Instantiate(
ps.get("input_type") + "InputHandler"));
Inputs.back().handler->Initialize(ps);
- return (Inputs.size() - 1);
+ return iid;
}
InputManager::Input_id_t
InputManager::GetInputId(std::string const &file_name) const {
for (size_t i = 0; i < Inputs.size(); ++i) {
if (Inputs[i].name == file_name) {
return i;
}
}
throw unknown_input() << "[ERROR]: Input file " << std::quoted(file_name)
<< " has not been loaded.";
}
IInputHandler const &
InputManager::GetInputHandler(InputManager::Input_id_t id) const {
if (id >= Inputs.size()) {
throw unknown_input() << "[ERROR]: Attempted to get input with id " << id
<< ", but only have " << Inputs.size()
<< " loaded inputs.";
}
return *Inputs[id].handler.get();
}
} // namespace input
} // namespace nuis
diff --git a/src/parameters/ParameterManager.cxx b/src/parameters/ParameterManager.cxx
index 76b17d6..4f39e47 100644
--- a/src/parameters/ParameterManager.cxx
+++ b/src/parameters/ParameterManager.cxx
@@ -1,133 +1,146 @@
#include "parameters/ParameterManager.hxx"
+#include "fhiclcpp/ParameterSet.h"
+
namespace nuis {
namespace params {
+ParameterManager::NamedParameter::NamedParameter()
+ : name(""), type(""), value(kDefaultValue), start(kDefaultValue),
+ min(kDefaultValue), max(kDefaultValue), step(kDefaultValue),
+ Penalty([](double) -> double { return 0; }) {}
+
+ParameterManager::NamedParameter::NamedParameter(NamedParameter &&other)
+ : name(std::move(other.name)), type(std::move(other.type)),
+ value(other.value), start(other.start), min(other.min), max(other.max),
+ step(other.step), Penalty(std::move(other.Penalty)) {}
+
ParameterManager *ParameterManager::_global_inst = nullptr;
ParameterManager::ParameterManager() : locked(false) {}
ParameterManager &ParameterManager::Get() {
if (!_global_inst) {
_global_inst = new ParameterManager();
}
return *_global_inst;
}
void ParameterManager::ValidateParamId(paramId_t pid) {
if (pid >= Parameters.size()) {
throw invalid_parameter_id()
<< "[ERROR]: Passed parameter id " << pid
<< ", but the ParameterManager only knows about " << Parameters.size()
<< " parameters.";
}
}
void ParameterManager::LockParameterList() { locked = true; }
void ParameterManager::UnlockParameterList() { locked = false; }
paramId_t
ParameterManager::EnsureParameterRegistered(fhicl::ParameterSet const &ps) {
if (locked) {
throw parameter_list_is_locked()
<< "[ERROR]: Attempted to register parameter: " << ps.to_string()
<< " when global ParameterManager was locked in state: "
<< StateString();
}
NamedParameter np;
np.name = ps.get("name");
np.type = ps.get("type");
paramId_t pid = GetParameterId(np.name, np.type);
if (pid != kParamUnhandled) {
return pid;
}
np.start = ps.get("start");
np.value = np.start;
np.min = ps.get("min", kDefaultLimit);
np.max = ps.get("max", kDefaultLimit);
np.step = ps.get("step");
pid = Parameters.size();
- Parameters.push_back(np);
+ Parameters.emplace_back(std::move(np));
return pid;
}
paramId_t ParameterManager::GetParameterId(std::string const &name,
std::string const &type) {
paramId_t pid = kParamUnhandled;
for (size_t p_it = 0; p_it < Parameters.size(); ++p_it) {
if (name != Parameters[p_it].name) {
continue;
}
if (type.size() && (type != Parameters[p_it].type)) {
continue;
}
// matches search, check if it is the first to match the search.
if (pid != kParamUnhandled) {
throw ambiguous_parameter_specified()
<< "[ERROR]: When searching for parameter by name-only, found at "
"least two matching parameters: { PID: "
<< pid << ", name: " << Parameters[pid].name
<< ", type: " << Parameters[pid].type << " } and { PID: " << p_it
<< ", name: " << Parameters[p_it].name
<< ", type: " << Parameters[p_it].type << " }";
}
pid = p_it;
}
return pid;
}
void ParameterManager::SetParameterValue(paramId_t pid, double val) {
ValidateParamId(pid);
if (!IsValidParameterValue(pid, val)) {
throw param_value_out_of_bounds()
- << "[ERROR]: Attempting to set parameter { PID: " << pid ", name: "
- << Parameters[pid].name << ", type: " << Parameters[pid].type
- << " } to " << val << ", but this is out of the allowed range ["
+ << "[ERROR]: Attempting to set parameter { PID: " << pid
+ << ", name: " << Parameters[pid].name
+ << ", type: " << Parameters[pid].type << " } to " << val
+ << ", but this is out of the allowed range ["
<< ((Parameters[pid].min == kDefaultLimit)
? "unbounded"
: std::to_string(Parameters[pid].min))
<< ","
<< ((Parameters[pid].max == kDefaultLimit)
? "unbounded"
: std::to_string(Parameters[pid].max))
<< "]";
}
Parameters[pid].value = val;
}
double ParameterManager::GetParameterValue(paramId_t pid) {
ValidateParamId(pid);
return Parameters[pid].value;
}
double ParameterManager::GetParameterStep(paramId_t pid) {
ValidateParamId(pid);
return Parameters[pid].step;
}
double ParameterManager::GetParameterStart(paramId_t pid) {
ValidateParamId(pid);
return Parameters[pid].start;
}
double ParameterManager::GetParameterMin(paramId_t pid) {
ValidateParamId(pid);
return Parameters[pid].min;
}
double ParameterManager::GetParameterMax(paramId_t pid) {
ValidateParamId(pid);
return Parameters[pid].max;
}
bool ParameterManager::IsValidParameterValue(paramId_t pid, double val) {
ValidateParamId(pid);
return (
((Parameters[pid].min == kDefaultLimit) || (Parameters[pid].min < val)) &&
((Parameters[pid].max == kDefaultLimit) || (Parameters[pid].max > val)));
}
std::string ParameterManager::StateString() {
return "Parameter Manager state:";
}
} // namespace params
} // namespace nuis
diff --git a/src/parameters/ParameterManager.hxx b/src/parameters/ParameterManager.hxx
index 07f8651..ae1de88 100644
--- a/src/parameters/ParameterManager.hxx
+++ b/src/parameters/ParameterManager.hxx
@@ -1,91 +1,97 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef PARAMETERS_PARAMETERMANAGER_HXX_SEEN
#define PARAMETERS_PARAMETERMANAGER_HXX_SEEN
#include "exception/exception.hxx"
#include
#include
#include
+#include
namespace fhicl {
class ParameterSet;
}
namespace nuis {
namespace params {
typedef size_t paramId_t;
static paramId_t const kParamUnhandled = std::numeric_limits::max();
static double const kDefaultLimit = 0xdeadbeef;
+static double const kDefaultValue = 0xdeadbeef;
class ParameterManager {
struct NamedParameter {
+ NamedParameter();
+ NamedParameter(NamedParameter const &) = delete;
+ NamedParameter(NamedParameter &&);
std::string name;
std::string type;
double value;
double start;
double min;
double max;
double step;
+ std::function Penalty;
};
std::vector Parameters;
// TMatrixD describing parameter covariance.
ParameterManager();
static ParameterManager *_global_inst;
void ValidateParamId(paramId_t);
bool locked;
public:
static ParameterManager &Get();
NEW_NUIS_EXCEPT(invalid_parameter_id);
NEW_NUIS_EXCEPT(param_value_out_of_bounds);
NEW_NUIS_EXCEPT(ambiguous_parameter_specified);
NEW_NUIS_EXCEPT(parameter_list_is_locked);
///\brief Lock the parameter list so that subsequent calls to
///EnsureParameterRegistered cause an exception to be thrown.
///
/// Useful for ensuring that plugins do not attempt to add parameters mid-fit.
void LockParameterList();
void UnlockParameterList();
paramId_t EnsureParameterRegistered(fhicl::ParameterSet const &);
paramId_t GetParameterId(std::string const &, std::string const &type = "");
void SetParameterValue(paramId_t, double);
double GetParameterValue(paramId_t);
double GetParameterStep(paramId_t);
double GetParameterStart(paramId_t);
double GetParameterMin(paramId_t);
double GetParameterMax(paramId_t);
bool IsValidParameterValue(paramId_t, double);
std::string StateString();
};
} // namespace params
} // namespace nuis
#endif
diff --git a/src/persistency/ROOTOutput.cxx b/src/persistency/ROOTOutput.cxx
index 7d83df5..8103ee8 100644
--- a/src/persistency/ROOTOutput.cxx
+++ b/src/persistency/ROOTOutput.cxx
@@ -1,45 +1,64 @@
#include "persistency/ROOTOutput.hxx"
#include "utility/ROOTUtility.hxx"
#include "config/GlobalConfiguration.hxx"
#include "fhiclcpp/ParameterSet.h"
#include "TFile.h"
namespace nuis {
namespace persistency {
struct NamedTFile {
std::string name;
- std::shared_ptr file;
+ std::unique_ptr file;
+
+ NamedTFile() : name(""), file(nullptr) {}
+ NamedTFile(NamedTFile &&other)
+ : name(std::move(other.name)), file(std::move(other.file)) {}
+ ~NamedTFile() {
+ if (file) {
+ file->Write();
+ file->Close();
+ }
+ }
};
+
static std::vector Files;
-std::shared_ptr GetOutputFile(std::string const &name) {
+std::unique_ptr &GetOutputFile(std::string const &name) {
for (NamedTFile &file : Files) {
if (file.name == name) {
return file.file;
}
}
fhicl::ParameterSet const &persistency =
config::GetDocument().get("persistency");
std::string file_name = persistency.get(name + ".output_file");
std::string open_opts =
persistency.get(name + ".open_mode", "CREATE");
NamedTFile ntf;
ntf.name = name;
- ntf.file = std::make_shared(file_name.c_str(), open_opts.c_str());
+ ntf.file = std::make_unique(file_name.c_str(), open_opts.c_str());
if (!ntf.file || !ntf.file->IsOpen()) {
throw utility::failed_to_open_TFile()
<< "[ERROR]: Failed to open output file: " << std::quoted(file_name)
<< " in write mode with opts = " << std::quoted(open_opts);
}
- Files.push_back(ntf);
+ Files.push_back(std::move(ntf));
return Files.back().file;
}
+
+void CloseOpenTFiles() {
+ for (NamedTFile &f : Files) {
+ std::cout << "[INFO]: Closing open TFile: " << f.name << " "
+ << f.file->GetName() << std::endl;
+ }
+ Files.clear();
+}
} // namespace persistency
} // namespace nuis
diff --git a/src/persistency/ROOTOutput.hxx b/src/persistency/ROOTOutput.hxx
index 3b45a71..1aa5bb6 100644
--- a/src/persistency/ROOTOutput.hxx
+++ b/src/persistency/ROOTOutput.hxx
@@ -1,66 +1,76 @@
#ifndef PERSITENCY_ROOTOUTPUT_HXX_SEEN
#define PERSITENCY_ROOTOUTPUT_HXX_SEEN
#include "exception/exception.hxx"
#include "TFile.h"
#include
#include
#include
namespace nuis {
namespace persistency {
NEW_NUIS_EXCEPT(WriteToOutputFile_nullptr);
/// Will get/open a TFile that is described in the global config
///
/// The named streams will be used to configure the file name and open mode from
/// the global config element persistency.: {file: output.root opts:
/// CREATE}
-std::shared_ptr GetOutputFile(std::string const &name = "default");
+std::unique_ptr &GetOutputFile(std::string const &name = "default");
template
inline void WriteToOutputFile(T *object, std::string const &object_name,
- std::string dir_name = "",
- std::string const &file_name = "default") {
+ std::string dir_name = "",
+ std::string const &file_name = "default") {
-if(!object){
- throw WriteToOutputFile_nullptr();
-}
+ if (!object) {
+ throw WriteToOutputFile_nullptr();
+ }
TDirectory *ogdir = gDirectory;
- std::shared_ptr f = GetOutputFile(file_name);
+ std::unique_ptr &f = GetOutputFile(file_name);
TDirectory *d = f.get();
while (dir_name.length()) {
size_t next_slash = dir_name.find_first_of('/');
std::string next_dir = dir_name.substr(0, next_slash);
if (next_slash != std::string::npos) {
dir_name = dir_name.substr(next_slash + 1);
} else {
dir_name = "";
}
TDirectory *nd = d->GetDirectory(next_dir.c_str());
if (!nd) {
nd = d->mkdir(next_dir.c_str());
}
nd->cd();
d = nd;
}
d->cd();
object->Write(object_name.c_str(), TObject::kOverwrite);
if (ogdir) {
ogdir->cd();
}
}
+
+template
+inline void WriteToOutputFile(std::unique_ptr &object,
+ std::string const &object_name,
+ std::string dir_name = "",
+ std::string const &file_name = "default") {
+ return WriteToOutputFile(object.get(), object_name, dir_name, file_name);
+}
+
+void CloseOpenTFiles();
} // namespace persistency
} // namespace nuis
#endif
diff --git a/src/plugins/CMakeLists.txt b/src/plugins/CMakeLists.txt
index 5c98f59..b2ab67f 100644
--- a/src/plugins/CMakeLists.txt
+++ b/src/plugins/CMakeLists.txt
@@ -1,5 +1,13 @@
+SET(plugins_implementation_files
+ LoadedPlugins.cxx)
+
set(plugins_header_files
-Instantiate.hxx
-traits.hxx)
+ Instantiate.hxx
+ NamedSO.hxx
+ traits.hxx)
+
+add_library(nuis_plugins SHARED ${plugins_implementation_files})
+target_link_libraries(nuis_plugins)
+install(TARGETS nuis_plugins DESTINATION lib)
install(FILES ${plugins_header_files} DESTINATION include/plugins)
diff --git a/src/plugins/Instantiate.hxx b/src/plugins/Instantiate.hxx
index 05c2bbd..8e7bf57 100644
--- a/src/plugins/Instantiate.hxx
+++ b/src/plugins/Instantiate.hxx
@@ -1,273 +1,254 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef PLUGINS_PLUGINMANAGER_HXX_SEEN
#define PLUGINS_PLUGINMANAGER_HXX_SEEN
+#include "plugins/NamedSO.hxx"
#include "plugins/traits.hxx"
#include "config/GlobalConfiguration.hxx"
#include "utility/FileSystemUtility.hxx"
#include "exception/exception.hxx"
#include "fhiclcpp/ParameterSet.h"
#include "string_parsers/to_string.hxx"
// linux
#include
#include
#include
#include
#include
#include
#include
// #define DEBUG_INSTANTIATE
namespace nuis {
namespace plugins {
+extern std::vector LoadedSharedObjects;
+
NEW_NUIS_EXCEPT(failed_to_find_instantiator);
NEW_NUIS_EXCEPT(malformed_plugin_interface);
NEW_NUIS_EXCEPT(failed_to_load_so);
typedef void *(*inst_fcn)();
typedef void (*dltr_fcn)(void *);
template struct PluginInstantiator {
std::string FQ_so_path;
std::string Base_classname;
std::string Classname;
void *dllib;
inst_fcn Instantiator;
dltr_fcn Deleter;
PluginInstantiator()
: FQ_so_path(""), Base_classname(""), Classname(""), dllib(nullptr),
Instantiator(nullptr), Deleter(nullptr) {}
PluginInstantiator(PluginInstantiator const &) = delete;
PluginInstantiator(PluginInstantiator &&other) {
FQ_so_path = std::move(other.FQ_so_path);
Base_classname = std::move(other.Base_classname);
Classname = std::move(other.Classname);
dllib = other.dllib;
Instantiator = other.Instantiator;
Deleter = other.Deleter;
other.FQ_so_path = "";
other.Base_classname = "";
other.Classname = "";
other.dllib = nullptr;
other.Instantiator = nullptr;
other.Deleter = nullptr;
}
typename plugin_traits::unique_ptr_t Instantiate() {
T *inst = reinterpret_cast((*Instantiator)());
dltr_fcn dltr = Deleter;
std::string cln = Classname;
std::function deleter = [=](T *inst) {
#ifdef DEBUG_INSTANTIATE
std::cout << "[INFO]: Deleting instance of " << cln << " with "
<< (void *)dltr << std::endl;
#endif
(*dltr)(inst);
};
return typename plugin_traits::unique_ptr_t(inst, deleter);
}
};
-struct NamedSO {
- std::string name;
- void *dllib;
-
- NamedSO() : name(""), dllib(nullptr) {}
- NamedSO(NamedSO const &) = delete;
- NamedSO(NamedSO &&other) : name(std::move(other.name)), dllib(other.dllib) {
- other.dllib = nullptr;
- }
-
- ~NamedSO() {
- if (dllib) {
-#ifdef DEBUG_INSTANTIATE
- std::cout << "[INFO]: dlclose on shared object: " << std::quoted(name)
- << std::endl;
-#endif
- dlclose(dllib);
- }
- }
-};
-
NamedSO &GetSharedObject(std::string const &FQPath) {
- static std::vector LoadedSharedObjects;
for (NamedSO &so : LoadedSharedObjects) {
if (so.name == FQPath) {
return so;
}
}
NamedSO so;
so.name = FQPath;
so.dllib = dlopen(FQPath.c_str(), RTLD_NOW | RTLD_GLOBAL);
char const *dlerr_cstr = dlerror();
std::string dlerr;
if (dlerr_cstr) {
dlerr = dlerr_cstr;
}
if (dlerr.length()) {
throw failed_to_load_so()
<< "[INFO]: Failed to load shared object: " << FQPath
<< " with dlerror: " << dlerr;
} else {
#ifdef DEBUG_INSTANTIATE
std::cout << "[INFO]: Loaded shared object " << FQPath << std::endl;
#endif
}
LoadedSharedObjects.push_back(std::move(so));
return LoadedSharedObjects.back();
}
template
typename plugin_traits::unique_ptr_t
Instantiate(std::string const &classname) {
static std::vector> LoadedPlugins;
fhicl::ParameterSet const &plugins =
config::GetDocument().get("plugins");
fhicl::ParameterSet const &search_paths =
plugins.get("search_paths");
std::vector plugin_search_dirs;
// Look for plugin search paths in sequence elements of the
// plugins.search_paths table
for (std::string const &key : search_paths.get_names()) {
if (!search_paths.is_key_to_sequence(key)) {
continue;
}
for (std::string const &path :
search_paths.get>(key)) {
plugin_search_dirs.push_back(path);
}
}
for (std::string path : plugin_search_dirs) {
path = utility::EnsureTrailingSlash(path);
for (std::string const &so_name :
utility::GetMatchingFiles(path, ".*\\.so")) {
for (PluginInstantiator &plugin : LoadedPlugins) {
if (plugin.FQ_so_path == (path + so_name) &&
(plugin.Base_classname == plugin_traits::interface_name()) &&
(plugin.Classname == classname)) {
#ifdef DEBUG_INSTANTIATE
std::cout << "[INFO]: Using already loaded PluginInstantiator"
<< std::endl;
#endif
return plugin.Instantiate();
}
}
PluginInstantiator plugin;
plugin.FQ_so_path = path + so_name;
plugin.Base_classname = plugin_traits::interface_name();
plugin.Classname = classname;
plugin.dllib = GetSharedObject(plugin.FQ_so_path).dllib;
char const *dlerr_cstr = nullptr;
std::string dlerr("");
plugin.Instantiator = reinterpret_cast(dlsym(
plugin.dllib,
plugin_traits::instantiator_function_name(classname).c_str()));
dlerr_cstr = dlerror();
if (dlerr_cstr) {
dlerr = dlerr_cstr;
}
if (dlerr_cstr) {
#ifdef DEBUG_INSTANTIATE
std::cout << "[INFO]: Failed to load appropriate instantiator method: "
<< plugin_traits::instantiator_function_name(classname)
<< " from shared object " << plugin.FQ_so_path;
#endif
continue;
} else {
#ifdef DEBUG_INSTANTIATE
std::cout << "[INFO]: Loaded instantiator method: "
<< plugin_traits::instantiator_function_name(classname)
<< " from shared object " << plugin.FQ_so_path << std::endl;
#endif
}
plugin.Deleter = reinterpret_cast(
dlsym(plugin.dllib,
plugin_traits::deleter_function_name(classname).c_str()));
dlerr_cstr = dlerror();
if (dlerr_cstr) {
dlerr = dlerr_cstr;
}
if (dlerr_cstr) {
throw malformed_plugin_interface()
<< "[ERROR]: Failed to load appropriate deleter method: "
<< plugin_traits::deleter_function_name(classname)
<< " from shared object " << plugin.FQ_so_path
<< " with error: " << std::quoted(dlerr);
} else {
#ifdef DEBUG_INSTANTIATE
std::cout << "[INFO]: Loaded deleter method: "
<< plugin_traits::deleter_function_name(classname)
<< " from shared object " << plugin.FQ_so_path << std::endl;
#endif
}
#ifdef DEBUG_INSTANTIATE
std::cout << "[INFO]: Checking if shared object "
<< std::quoted(plugin.FQ_so_path)
<< " knows how to instantiate class " << std::quoted(classname)
<< " via interface "
<< std::quoted(plugin_traits::interface_name()) << std::endl;
#endif
LoadedPlugins.push_back(std::move(plugin));
return LoadedPlugins.back().Instantiate();
}
}
throw failed_to_find_instantiator()
<< "[ERROR]: Failed to find instantiator for classname: "
<< std::quoted(classname) << " using interface "
<< std::quoted(plugin_traits::interface_name())
<< " from configured search paths: "
<< fhicl::string_parsers::T2Str>(
plugin_search_dirs);
}
} // namespace plugins
} // namespace nuis
#endif
diff --git a/src/plugins/LoadedPlugins.cxx b/src/plugins/LoadedPlugins.cxx
new file mode 100644
index 0000000..86e1edf
--- /dev/null
+++ b/src/plugins/LoadedPlugins.cxx
@@ -0,0 +1,9 @@
+#include "plugins/NamedSO.hxx"
+
+#include
+
+namespace nuis {
+namespace plugins {
+std::vector LoadedSharedObjects;
+}
+} // namespace nuis
diff --git a/src/plugins/NamedSO.hxx b/src/plugins/NamedSO.hxx
new file mode 100644
index 0000000..5137700
--- /dev/null
+++ b/src/plugins/NamedSO.hxx
@@ -0,0 +1,32 @@
+#ifndef PLUGINS_NAMEDSO_HXX_SEEN
+#define PLUGINS_NAMEDSO_HXX_SEEN
+
+// linux
+#include
+
+#include
+#include
+#include
+
+struct NamedSO {
+ std::string name;
+ void *dllib;
+
+ NamedSO() : name(""), dllib(nullptr) {}
+ NamedSO(NamedSO const &) = delete;
+ NamedSO(NamedSO &&other) : name(std::move(other.name)), dllib(other.dllib) {
+ other.dllib = nullptr;
+ }
+
+ ~NamedSO() {
+ if (dllib) {
+#ifdef DEBUG_INSTANTIATE
+ std::cout << "[INFO]: dlclose on shared object: " << std::quoted(name)
+ << std::endl;
+#endif
+ dlclose(dllib);
+ }
+ }
+};
+
+#endif
diff --git a/src/samples/CMakeLists.txt b/src/samples/CMakeLists.txt
index ed8bfaf..a31710c 100644
--- a/src/samples/CMakeLists.txt
+++ b/src/samples/CMakeLists.txt
@@ -1,20 +1,23 @@
set(samples_header_files
- ISample.hxx
+ IEventProcessor.hxx
IDataComparison.hxx
SimpleDataComparison.hxx)
install(FILES ${samples_header_files} DESTINATION include/samples)
add_subdirectory(MCTools)
add_subdirectory(nuA)
cmessage(DEBUG "INuADataComparisons: ${INuADataComparisons}")
SET(INuADataComparisons_List)
if(NOT IDataComparisons STREQUAL "")
- string(REPLACE ";" "," INuADataComparisons_List ${INuADataComparisons})
+ string(REPLACE ";" ", " INuADataComparisons_List "${INuADataComparisons}")
endif()
+cmessage(DEBUG "INuADataComparisons_List: ${INuADataComparisons_List}")
+
+
SET(INuADataComparisons_List ${INuADataComparisons_List} PARENT_SCOPE)
SET(INuADataComparisons_FHiCL ${INuADataComparisons_FHiCL} PARENT_SCOPE)
diff --git a/src/samples/IDataComparison.hxx b/src/samples/IDataComparison.hxx
index cfda1fc..1fdc67d 100644
--- a/src/samples/IDataComparison.hxx
+++ b/src/samples/IDataComparison.hxx
@@ -1,73 +1,73 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef SAMPLES_IDATACOMPARISON_HXX_SEEN
#define SAMPLES_IDATACOMPARISON_HXX_SEEN
-#include "samples/ISample.hxx"
+#include "samples/IEventProcessor.hxx"
#include "fhiclcpp/ParameterSet.h"
#include
#include
-class IDataComparison : public ISample {
+class IDataComparison : public IEventProcessor {
public:
virtual double GetGOF() = 0;
virtual double GetNDOGuess() = 0;
virtual std::string GetJournalReference() {
std::stringstream ss("");
ss << "Unknown Journal Ref. for IDataComparison: " << std::quoted(Name());
return ss.str();
}
virtual std::string GetTargetMaterial() {
std::stringstream ss("");
ss << "Unknown Target material for IDataComparison: " << std::quoted(Name());
return ss.str();
}
virtual std::string GetFluxDescription() {
std::stringstream ss("");
ss << "Unknown Flux description for IDataComparison: " << std::quoted(Name());
return ss.str();
}
virtual std::string GetSignalDescription() {
std::stringstream ss("");
ss << "Unknown Signal description for IDataComparison: " << std::quoted(
Name());
return ss.str();
}
virtual std::string GetDocumentation() {
std::stringstream ss("");
ss << "No documentation provided for IDataComparison: " << std::quoted(Name());
return ss.str();
}
virtual fhicl::ParameterSet GetExampleConfiguration() {
return fhicl::ParameterSet();
}
};
DECLARE_PLUGIN_INTERFACE(IDataComparison);
#endif
diff --git a/src/samples/ISample.hxx b/src/samples/IEventProcessor.hxx
similarity index 74%
rename from src/samples/ISample.hxx
rename to src/samples/IEventProcessor.hxx
index 0b195b6..90eff1a 100644
--- a/src/samples/ISample.hxx
+++ b/src/samples/IEventProcessor.hxx
@@ -1,113 +1,119 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
-#ifndef SAMPLES_ISAMPLE_HXX_SEEN
-#define SAMPLES_ISAMPLE_HXX_SEEN
+#ifndef SAMPLES_IEventProcessor_HXX_SEEN
+#define SAMPLES_IEventProcessor_HXX_SEEN
#include "plugins/traits.hxx"
#include "exception/exception.hxx"
#include "config/GlobalConfiguration.hxx"
#include "fhiclcpp/ParameterSet.h"
+#include "TH1.h"
+
#include
#include
#include
namespace nuis {
namespace event {
class FullEvent;
class MinimalEvent;
} // namespace event
} // namespace nuis
-#define ISAMPLE_DEBUG(v) \
+#define IEventProcessor_DEBUG(v) \
if (verb > 2) { \
std::cout << "[DEBUG]: " << v << std::endl; \
}
-#define ISAMPLE_INFO(v) \
+#define IEventProcessor_INFO(v) \
if (verb > 1) { \
std::cout << "[INFO]: " << v << std::endl; \
}
-#define ISAMPLE_WARN(v) \
+#define IEventProcessor_WARN(v) \
if (verb > 0) { \
std::cout << "[WARN] " << __FILENAME__ << ":" << __LINE__ << " : " << v \
<< std::endl; \
}
-class ISample {
+class IEventProcessor {
protected:
- NEW_NUIS_EXCEPT(unknown_ISample_verbosity);
+ NEW_NUIS_EXCEPT(unknown_IEventProcessor_verbosity);
enum sample_verbosity { kSilent = 0, kWarn = 1, kReticent = 2, kVerbose = 3 };
sample_verbosity verb;
void SetSampleVerbosity(std::string v) {
v = nuis::config::GetDocument().get(
"global.sample.verbosity_override", v);
if (v == "Silent") {
verb = kSilent;
} else if (v == "Warn") {
verb = kWarn;
} else if (v == "Reticent") {
verb = kReticent;
} else if (v == "Verbose") {
verb = kVerbose;
} else {
- throw unknown_ISample_verbosity()
- << "[ERROR]: Failed to parse ISample verbosity setting from: "
+ throw unknown_IEventProcessor_verbosity()
+ << "[ERROR]: Failed to parse IEventProcessor verbosity setting from: "
<< std::quoted(v);
}
}
public:
- NEW_NUIS_EXCEPT(uninitialized_ISample);
- NEW_NUIS_EXCEPT(unimplemented_ISample_optional_method);
+ NEW_NUIS_EXCEPT(uninitialized_IEventProcessor);
+ NEW_NUIS_EXCEPT(unimplemented_IEventProcessor_optional_method);
- ISample() {
+ IEventProcessor() {
SetSampleVerbosity(nuis::config::GetDocument().get(
"global.sample.verbosity_default", "Silent"));
+
+ TH1::SetDefaultSumw2();
}
virtual void Initialize(fhicl::ParameterSet const &) = 0;
// Interface for processing a single, external event
//
- // ISamples are not required to implement processing events from 'outside'.
+ // IEventProcessors are not required to implement processing events from
+ // 'outside'.
virtual void ProcessEvent(nuis::event::FullEvent const &) {
- throw unimplemented_ISample_optional_method()
- << "[ERROR]: ISample " << Name() << " does not implement ProcessEvent.";
+ throw unimplemented_IEventProcessor_optional_method()
+ << "[ERROR]: IEventProcessor " << Name()
+ << " does not implement ProcessEvent.";
}
virtual void
ProcessSample(size_t nmax = std::numeric_limits::max()) = 0;
virtual void Write() = 0;
virtual std::string Name() = 0;
- virtual ~ISample() {}
+ virtual ~IEventProcessor() {}
};
-DECLARE_PLUGIN_INTERFACE(ISample);
+DECLARE_PLUGIN_INTERFACE(IEventProcessor);
#endif
diff --git a/src/samples/MCTools/NuisToNuWro.cxx b/src/samples/MCTools/NuisToNuWro.cxx
index b408c53..dd1e6c4 100644
--- a/src/samples/MCTools/NuisToNuWro.cxx
+++ b/src/samples/MCTools/NuisToNuWro.cxx
@@ -1,103 +1,103 @@
-#include "samples/ISample.hxx"
+#include "samples/IEventProcessor.hxx"
#include "event/FullEvent.hxx"
#include "input/InputManager.hxx"
#include "utility/FullEventUtility.hxx"
#include "utility/ROOTUtility.hxx"
#include "fhiclcpp/ParameterSet.h"
#include "generator/utility/NuWroUtility.hxx"
#include
#include
using namespace nuis::event;
using namespace nuis::input;
using namespace nuis::utility;
using namespace nuis::nuwrotools;
-class NuisToNuWro : public ISample {
+class NuisToNuWro : public IEventProcessor {
public:
InputManager::Input_id_t fIH_id;
- std::unique_ptr fOutputTree;
+ TreeFile fOutputTree;
NuWroEvent *fOutputEvent;
NuisToNuWro()
: fIH_id(std::numeric_limits::max()),
- fOutputTree(nullptr), fOutputEvent(nullptr) {}
+ fOutputTree(), fOutputEvent(nullptr) {}
void Initialize(fhicl::ParameterSet const &ps) {
fIH_id = InputManager::Get().EnsureInputLoaded(ps);
fOutputTree =
MakeNewTTree(ps.get("output_file"), "treeout", "RECREATE");
- fOutputTree->tree->Branch("e", &fOutputEvent);
+ fOutputTree.tree->Branch("e", &fOutputEvent);
}
void ProcessEvent(FullEvent const &ps) {
fOutputEvent->in.clear();
fOutputEvent->out.clear();
fOutputEvent->post.clear();
std::pair NuMode = GetFlagsDynEquivalent(ps.mode);
fOutputEvent->flag = NuMode.first;
fOutputEvent->dyn = NuMode.second;
for (Particle const &part : GetISParticles(ps)) {
particle nuwro_part(part.pdg, part.P4.M());
nuwro_part.set_momentum(vect(part.P4.X(), part.P4.Y(), part.P4.Z()));
fOutputEvent->in.push_back(nuwro_part);
}
for (Particle const &part : GetPrimaryFSParticles(ps)) {
particle nuwro_part(part.pdg, part.P4.M());
nuwro_part.set_momentum(vect(part.P4.X(), part.P4.Y(), part.P4.Z()));
fOutputEvent->out.push_back(nuwro_part);
}
for (Particle const &part : GetNuclearLeavingParticles(ps)) {
particle nuwro_part(part.pdg, part.P4.M());
nuwro_part.set_momentum(vect(part.P4.X(), part.P4.Y(), part.P4.Z()));
fOutputEvent->post.push_back(nuwro_part);
}
- fOutputTree->tree->Fill();
+ fOutputTree.tree->Fill();
}
void ProcessSample(size_t nmax) {
if (fIH_id == std::numeric_limits::max()) {
- throw uninitialized_ISample();
+ throw uninitialized_IEventProcessor();
}
IInputHandler const &IH = InputManager::Get().GetInputHandler(fIH_id);
size_t NEvsToProcess = std::min(nmax, IH.GetNEvents());
size_t NToShout = NEvsToProcess / 10;
std::cout << "[INFO]: Processing " << NEvsToProcess
<< " input events to NuWro format." << std::endl;
size_t n = 0;
for (auto const &fe : IH) {
if (++n > NEvsToProcess) {
break;
}
if (NToShout && !(n % NToShout)) {
std::cout << "[INFO]: Processed " << n << "/" << NEvsToProcess
<< " events." << std::endl;
}
ProcessEvent(fe);
}
}
- void Write() { fOutputTree->file->Write(); }
+ void Write() { fOutputTree.file->Write(); }
std::string Name() { return "NuisToNuWro"; }
};
-DECLARE_PLUGIN(ISample, NuisToNuWro);
+DECLARE_PLUGIN(IEventProcessor, NuisToNuWro);
diff --git a/src/samples/MCTools/VerboseEventSummary.cxx b/src/samples/MCTools/VerboseEventSummary.cxx
index 726a5a5..a007433 100644
--- a/src/samples/MCTools/VerboseEventSummary.cxx
+++ b/src/samples/MCTools/VerboseEventSummary.cxx
@@ -1,56 +1,56 @@
-#include "samples/ISample.hxx"
+#include "samples/IEventProcessor.hxx"
#include "event/FullEvent.hxx"
#include "input/InputManager.hxx"
#include
#include
using namespace nuis::event;
using namespace nuis::input;
-class VerboseEventSummary : public ISample {
+class VerboseEventSummary : public IEventProcessor {
public:
InputManager::Input_id_t fIH_id;
VerboseEventSummary()
: fIH_id(std::numeric_limits::max()) {}
void Initialize(fhicl::ParameterSet const &ps) {
fIH_id = InputManager::Get().EnsureInputLoaded(ps);
}
void ProcessEvent(FullEvent const &ps) {
std::cout << ps.to_string() << std::endl;
}
void ProcessSample(size_t nmax) {
if (fIH_id == std::numeric_limits::max()) {
- throw uninitialized_ISample();
+ throw uninitialized_IEventProcessor();
}
IInputHandler const &IH = InputManager::Get().GetInputHandler(fIH_id);
size_t NEvsToProcess = std::min(nmax, IH.GetNEvents());
size_t NToShout = NEvsToProcess / 10;
std::cout << "[INFO]: Processing " << NEvsToProcess << " input events."
<< std::endl;
size_t n = 0;
for (auto const &fe : IH) {
if (++n > nmax) {
break;
}
if (NToShout && !(n % NToShout)) {
std::cout << "[INFO]: Processed " << n << "/" << NEvsToProcess
<< " events." << std::endl;
}
ProcessEvent(fe);
}
}
void Write() {}
std::string Name() { return "VerboseEventSummary"; }
};
-DECLARE_PLUGIN(ISample, VerboseEventSummary);
+DECLARE_PLUGIN(IEventProcessor, VerboseEventSummary);
diff --git a/src/samples/SimpleDataComparison.hxx b/src/samples/SimpleDataComparison.hxx
index bdce072..36f8e7a 100644
--- a/src/samples/SimpleDataComparison.hxx
+++ b/src/samples/SimpleDataComparison.hxx
@@ -1,375 +1,424 @@
// Copyright 2018 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
/*******************************************************************************
* This file is part of NUISANCE.
*
* NUISANCE is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* NUISANCE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with NUISANCE. If not, see .
*******************************************************************************/
#ifndef SAMPLES_SIMPLEDATACOMPARISON_HXX_SEEN
#define SAMPLES_SIMPLEDATACOMPARISON_HXX_SEEN
#include "samples/IDataComparison.hxx"
#include "event/FullEvent.hxx"
#include "input/InputManager.hxx"
#include "persistency/ROOTOutput.hxx"
#include "utility/FileSystemUtility.hxx"
#include "utility/HistogramUtility.hxx"
+#include "utility/KinematicUtility.hxx"
#include "utility/ROOTUtility.hxx"
#include "utility/StatsUtility.hxx"
#include
#include
#include
#include
#include
-template struct TH_dim_helper {};
-template <> struct TH_dim_helper<1> { typedef TH1D type; };
-template <> struct TH_dim_helper<2> { typedef TH2D type; };
-
-template class SimpleDataComparison : public IDataComparison {
+template ::type>
+class SimpleDataComparison : public IDataComparison {
NEW_NUIS_EXCEPT(invalid_SimpleDataComparison_initialization);
NEW_NUIS_EXCEPT(SimpleDataComparison_already_finalized);
protected:
+ using HistType = HT;
+ using TH_Help = typename nuis::utility::TH_Helper;
+
nuis::input::InputManager::Input_id_t fIH_id;
std::string write_directory;
size_t NMaxSample_override;
int fIsShapeOnly;
int fIsFluxUnfolded;
std::vector fSignalCache;
- std::vector> fProjectionCache;
+ std::vector> fProjectionCache;
std::string fDataInputDescriptor;
- std::unique_ptr::type> fData;
+ std::unique_ptr fData;
std::string fMaskInputDescriptor;
- std::unique_ptr::type> fMask;
+ std::unique_ptr fMask;
std::string fCovarianceInputDescriptor;
- std::unique_ptr fCovariance;
- std::unique_ptr::type> fPrediction;
- std::unique_ptr::type> fPrediction_xsec;
- std::unique_ptr::type> fPrediction_shape;
- std::unique_ptr::type> fPrediction_comparison;
+ std::unique_ptr fCovariance;
+ std::unique_ptr fPrediction;
+ std::unique_ptr fPrediction_xsec;
+ std::unique_ptr fPrediction_shape;
+ std::unique_ptr fPrediction_comparison;
bool fComparisonFinalized;
std::string fJournalReference;
std::string fTargetMaterial;
std::string fFluxDescription;
std::string fSignalDescription;
- std::pair EnuRange;
+ nuis::utility::ENuRange energy_cut;
std::function IsSigFunc;
- std::function(nuis::event::FullEvent const &)>
+ std::function(nuis::event::FullEvent const &)>
CompProjFunc;
+ // If assigned by subclass will be called on for all events, bool signifies
+ // whether the event was selected.
+ std::function
+ ProcessExtraFunc;
+
public:
SimpleDataComparison() {
fIH_id = std::numeric_limits::max();
write_directory = "";
NMaxSample_override = std::numeric_limits::max();
fDataInputDescriptor = "";
fData = nullptr;
fMaskInputDescriptor = "";
fMask = nullptr;
fCovarianceInputDescriptor = "";
fCovariance = nullptr;
fPrediction = nullptr;
fPrediction_xsec = nullptr;
fPrediction_shape = nullptr;
fPrediction_comparison = nullptr;
fComparisonFinalized = false;
IsSigFunc = [](nuis::event::FullEvent const &) -> bool { return true; };
CompProjFunc =
- [](nuis::event::FullEvent const &) -> std::array {
- std::array arr;
- for (double &el : arr) {
+ [](nuis::event::FullEvent const &) -> std::array {
+ std::array arr;
+ for (NumericT &el : arr) {
el = 0;
}
return arr;
};
+ ProcessExtraFunc =
+ std::function();
+
fJournalReference = "";
fTargetMaterial = "";
fFluxDescription = "";
fSignalDescription = "";
fIsShapeOnly = -1;
fIsFluxUnfolded = -1;
- EnuRange = std::pair{std::numeric_limits::max(),
+
+ energy_cut = nuis::utility::ENuRange{std::numeric_limits::max(),
std::numeric_limits::max()};
}
+ fhicl::ParameterSet fGlobalConfig;
+ fhicl::ParameterSet fInstanceConfig;
+
void ReadGlobalConfigDefaults() {
- fhicl::ParameterSet const &global_sample_configuration =
- nuis::config::GetDocument().get(
- std::string("global.sample_configuration.") + Name(),
- fhicl::ParameterSet());
+ fGlobalConfig = nuis::config::GetDocument().get(
+ std::string("global.sample_configuration.") + Name(),
+ fhicl::ParameterSet());
if (!fJournalReference.length()) {
- fJournalReference = global_sample_configuration.get(
- "journal_reference", fJournalReference);
+ fJournalReference = fGlobalConfig.get("journal_reference",
+ fJournalReference);
}
if (!fTargetMaterial.length()) {
- fTargetMaterial = global_sample_configuration.get(
- "target_material", fTargetMaterial);
+ fTargetMaterial =
+ fGlobalConfig.get("target_material", fTargetMaterial);
}
if (!fFluxDescription.length()) {
- fFluxDescription = global_sample_configuration.get(
- "flux_description", fFluxDescription);
+ fFluxDescription =
+ fGlobalConfig.get("flux_description", fFluxDescription);
}
if (!fSignalDescription.length()) {
- fSignalDescription = global_sample_configuration.get(
- "signal_description", fSignalDescription);
+ fSignalDescription = fGlobalConfig.get("signal_description",
+ fSignalDescription);
}
if (fIsShapeOnly == -1) {
- fIsShapeOnly = global_sample_configuration.get("shape_only", false);
+ fIsShapeOnly = fGlobalConfig.get("shape_only", false);
}
if (fIsFluxUnfolded == -1) {
- fIsFluxUnfolded =
- global_sample_configuration.get("flux_unfolded", false);
+ fIsFluxUnfolded = fGlobalConfig.get("flux_unfolded", false);
}
- if ((EnuRange.first == std::numeric_limits::max()) &&
- (global_sample_configuration.has_key("enu_range"))) {
- EnuRange = global_sample_configuration.get>(
- "enu_range");
+ if ((energy_cut.first == std::numeric_limits::max()) &&
+ (fGlobalConfig.has_key("enu_range"))) {
+ energy_cut = fGlobalConfig.get