diff --git a/.github/workflows/validation.yml b/.github/workflows/validation.yml index 1b02a2e0..9cc34f1e 100644 --- a/.github/workflows/validation.yml +++ b/.github/workflows/validation.yml @@ -4,14 +4,16 @@ on: workflow_dispatch: push: - branches: [master] + branches: [ master ] pull_request: - branches: [master] + branches: [ master ] + release: env: CMAKE_BUILD_TYPE: Release REST_PATH: /rest/framework/install REST_FRAMEWORK_SOURCE_DIR: /rest/framework/source + BRANCH_NAME: ${{ github.head_ref || github.ref_name }} defaults: run: @@ -22,23 +24,22 @@ jobs: name: Install framework with restG4 runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev steps: - name: Print version of dependencies for image run: version.sh - - name: Checkout framework run: | git clone https://github.com/rest-for-physics/framework.git ${{ env.REST_FRAMEWORK_SOURCE_DIR }} cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} - ./scripts/checkoutRemoteBranch.sh ${{ github.sha }} + ./scripts/checkoutRemoteBranch.sh ${{ env.BRANCH_NAME }} + git log -1 --stat - uses: actions/checkout@v3 - name: Setup, build and install run: | cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} git submodule init source/libraries/geant4 && git submodule update source/libraries/geant4 - cd source/libraries/geant4/ && ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/scripts/checkoutRemoteBranch.sh ${{ github.sha }} + cd source/libraries/geant4/ && ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/scripts/checkoutRemoteBranch.sh ${{ env.BRANCH_NAME }} cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} rm -rf source/packages/restG4/ && cp -r $GITHUB_WORKSPACE source/packages/restG4 cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} @@ -47,64 +48,40 @@ jobs: make -j4 install - name: Cache framework installation id: framework-install-restG4-cache - uses: actions/cache@v2 - with: - path: ${{ env.REST_PATH }} - key: ${{ github.sha }} - - framework-install-old-Geant4: - name: Install framework with restG4 on an older reference version of Geant4 - runs-on: ubuntu-latest - container: - image: ghcr.io/lobis/root-geant4-garfield:cpp17_ROOT-v6-26-00_Geant4-v10.4.3_Garfield-4ac1e62f - - steps: - - name: Print version of dependencies for image - run: version.sh - - - name: Checkout framework - run: | - git clone https://github.com/rest-for-physics/framework.git ${{ env.REST_FRAMEWORK_SOURCE_DIR }} - cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} - ./scripts/checkoutRemoteBranch.sh ${{ github.sha }} - - uses: actions/checkout@v3 - - name: Setup, build and install - run: | - cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} - git submodule init source/libraries/geant4 && git submodule update source/libraries/geant4 - cd source/libraries/geant4/ && ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/scripts/checkoutRemoteBranch.sh ${{ github.sha }} - cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} - rm -rf source/packages/restG4/ && cp -r $GITHUB_WORKSPACE source/packages/restG4 - cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} - mkdir -p ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/build && cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/build - cmake ../ -DCMAKE_BUILD_TYPE=${{ env.CMAKE_BUILD_TYPE }} -DREST_WELCOME=ON -DREST_G4=ON -DCMAKE_INSTALL_PREFIX=${{ env.REST_PATH }} - make -j4 install - - name: Cache framework installation - id: framework-install-restG4-cache-old-Geant4 - uses: actions/cache@v2 + uses: actions/cache@v3 with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }}-old-Geant4 + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} check-installation: name: Check framework and restG4 are accessible runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [framework-install] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ framework-install ] steps: - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} - + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} + - name: Check libraries + run: | + echo $LD_LIBRARY_PATH + ls -lht ${{ env.REST_PATH }} + ls -lht ${{ env.REST_PATH }}/lib | grep .so + ls -lht ${{ env.REST_PATH }}/bin + - name: Check root + run: | + root-config --version + root -b -q - name: Check framework run: | source ${{ env.REST_PATH }}/thisREST.sh - restRoot -b -q + cat ${{ env.REST_PATH }}/thisREST.sh + cat ${{ env.REST_PATH }}/bin/rest-config + rest-config --welcome - name: Check restG4 run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -114,19 +91,16 @@ jobs: name: Install restG4 as a standalone application using the installed framework runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [framework-install] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ framework-install, check-installation ] steps: - uses: actions/checkout@v3 - - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} - + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Build as standalone run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -136,19 +110,19 @@ jobs: /rest/restG4-standalone/install/bin/restG4 --help restG4-examples-01: + name: "Example 01: NLDBD" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -157,19 +131,19 @@ jobs: restRoot -b -q Validate.C'("Run00001_NLDBD_Test.root")' restG4-examples-03: + name: "Example 03: Fluorescence" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation, restG4-examples-01] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation, restG4-examples-01 ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -179,19 +153,19 @@ jobs: restRoot -b -q Validate.C'("gamma_fluorescence_analysis.root")' restG4-examples-04: + name: "Example 04: Cosmic Muons" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -204,19 +178,19 @@ jobs: restRoot -b -q ValidateCircle.C'("muons_from_circle.root")' restG4-examples-05: + name: "Example 05: PandaX-III" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation, restG4-examples-01] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation, restG4-examples-01 ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -225,19 +199,19 @@ jobs: restRoot -b -q ValidateG4.C'("Xe136bb0n_n2E06.root")' restG4-examples-06: + name: "Example 06: Ion recoils" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation, restG4-examples-01] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation, restG4-examples-01 ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -246,45 +220,45 @@ jobs: restRoot -b -q Validate.C'("Run00001_F20_Recoils.root")' restG4-examples-07: + name: "Example 07: Decay" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation, restG4-examples-01] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation, restG4-examples-01 ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh cd ${{ env.REST_PATH }}/examples/restG4/07.FullChainDecay/ restG4 fullChain.rml -o Run00001_U238_FullChainDecay.root - restRoot -b -q Validate.C'("Run00001_U238_FullChainDecay.root", 14)' + restRoot -b -q Validate.C'("Run00001_U238_FullChainDecay.root", 16)' restG4 singleDecay.rml -o Run00002_U238_SingleChainDecay.root restRoot -b -q Validate.C'("Run00002_U238_SingleChainDecay.root", 1)' - export REST_ISOTOPE="Be7" + export REST_ISOTOPE=Be7 restG4 singleDecay.rml restRoot -b -q Validate.C'("Run00002_Be7_SingleChainDecay.root", 1)' restG4-examples-08: + name: "Example 08: Alphas" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -306,19 +280,19 @@ jobs: restRoot -b -q Validate.C restG4-examples-09: + name: "Example 09: Shielding" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation, restG4-examples-01] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation, restG4-examples-01 ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -327,19 +301,19 @@ jobs: restRoot -b -q Validate.C'("Run00001_Pb210_Shielding.root")' restG4-examples-10: + name: "Example 10: Geometries" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation, restG4-examples-01] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation, restG4-examples-01 ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh @@ -348,23 +322,134 @@ jobs: restRoot -b -q Validate.C'("Run00001_Assembly_Assembly.root")' restG4-examples-11: + name: "Example 11: X-Rays" runs-on: ubuntu-latest container: - image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics - needs: [check-installation, restG4-examples-01] - + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev + needs: [ check-installation, restG4-examples-01 ] steps: - uses: actions/checkout@v3 - name: Restore cache - uses: actions/cache@v2 + uses: actions/cache@v3 id: framework-install-restG4-cache with: path: ${{ env.REST_PATH }} - key: ${{ github.sha }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} - name: Run example run: | source ${{ env.REST_PATH }}/thisREST.sh cd ${{ env.REST_PATH }}/examples/restG4/11.Xrays/ - restG4 xrays.rml -o RestG4_xrays_00001_ArgonISO.root - restManager --c analysis.rml --f RestG4_xrays_00001_ArgonISO.root - restRoot -b -q Validate.C'("Run_g4Analysis_00001_xrays.root")' + restG4 xrays.rml -o xrays_simulation.root + restManager --c analysis.rml --f xrays_simulation.root --o xrays_simulation_analysis.root + restRoot -b -q GetQE.C'("xrays_simulation_analysis.root")' + + # Reference version of Geant4 + + framework-install-reference: + name: Install framework with restG4 on an older reference version of Geant4 (v10.4.3) + runs-on: ubuntu-latest + container: + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-reference-jun2022 + steps: + - name: Print version of dependencies for image + run: version.sh + - name: Checkout framework + run: | + git clone https://github.com/rest-for-physics/framework.git ${{ env.REST_FRAMEWORK_SOURCE_DIR }} + cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} + ./scripts/checkoutRemoteBranch.sh ${{ env.BRANCH_NAME }} + git log -1 --stat + - uses: actions/checkout@v3 + - name: Setup, build and install + run: | + cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} + git submodule init source/libraries/geant4 && git submodule update source/libraries/geant4 + cd source/libraries/geant4/ && ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/scripts/checkoutRemoteBranch.sh ${{ env.BRANCH_NAME }} + cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} + rm -rf source/packages/restG4/ && cp -r $GITHUB_WORKSPACE source/packages/restG4 + cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }} + mkdir -p ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/build && cd ${{ env.REST_FRAMEWORK_SOURCE_DIR }}/build + cmake ../ -DCMAKE_BUILD_TYPE=${{ env.CMAKE_BUILD_TYPE }} -DREST_WELCOME=ON -DREST_G4=ON -DCMAKE_INSTALL_PREFIX=${{ env.REST_PATH }} + make -j4 install + - name: Cache framework installation + id: framework-install-restG4-cache-old-Geant4 + uses: actions/cache@v3 + with: + path: ${{ env.REST_PATH }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} + + check-installation-reference: + name: Check framework and restG4 are accessible + runs-on: ubuntu-latest + container: + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-reference-jun2022 + needs: [ framework-install-reference ] + steps: + - name: Restore cache + uses: actions/cache@v3 + id: framework-install-restG4-cache-reference + with: + path: ${{ env.REST_PATH }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} + - name: Check root + run: | + root-config --version + root -b -q + - name: Check framework + run: | + source ${{ env.REST_PATH }}/thisREST.sh + cat ${{ env.REST_PATH }}/thisREST.sh + cat ${{ env.REST_PATH }}/bin/rest-config + rest-config --welcome + - name: Check restG4 + run: | + source ${{ env.REST_PATH }}/thisREST.sh + restG4 --help + + + restG4-examples-01-reference: + name: "Reference Example 01: NLDBD" + runs-on: ubuntu-latest + container: + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-reference-jun2022 + needs: [ framework-install-reference, check-installation-reference ] + steps: + - uses: actions/checkout@v3 + - name: Restore cache + uses: actions/cache@v3 + id: framework-install-restG4-cache-reference + with: + path: ${{ env.REST_PATH }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} + - name: Run example + run: | + source ${{ env.REST_PATH }}/thisREST.sh + cd ${{ env.REST_PATH }}/examples/restG4/01.NLDBD/ + restG4 NLDBD.rml -o Run00001_NLDBD_Test.root + restRoot -b -q Validate.C'("Run00001_NLDBD_Test.root")' + + restG4-examples-07-reference: + name: "Reference Example 07: Decay" + runs-on: ubuntu-latest + container: + image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-reference-jun2022 + needs: [ framework-install-reference, check-installation-reference ] + steps: + - uses: actions/checkout@v3 + - name: Restore cache + uses: actions/cache@v3 + id: framework-install-restG4-cache-reference + with: + path: ${{ env.REST_PATH }} + key: ${{ env.BRANCH_NAME }}-${{ github.sha }} + - name: Run example + run: | + source ${{ env.REST_PATH }}/thisREST.sh + cd ${{ env.REST_PATH }}/examples/restG4/07.FullChainDecay/ + restG4 fullChain.rml -o Run00001_U238_FullChainDecay.root + restRoot -b -q Validate.C'("Run00001_U238_FullChainDecay.root", 15)' + restG4 singleDecay.rml -o Run00002_U238_SingleChainDecay.root + restRoot -b -q Validate.C'("Run00002_U238_SingleChainDecay.root", 1)' + export REST_ISOTOPE=Be7 + restG4 singleDecay.rml + restRoot -b -q Validate.C'("Run00002_Be7_SingleChainDecay.root", 1)' diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index aebb3cc8..b53b0796 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -1,10 +1,10 @@ -image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics +image: ghcr.io/lobis/root-geant4-garfield:rest-for-physics-dev stages: - - build - - build-standalone - - loadRESTLibs - - examples + - Build + - Build Standalone + - Load REST Libraries + - Run Examples before_script: - export USER="rest" @@ -18,7 +18,7 @@ before_script: # This is the one we use to run validations tests. We save the installation as an artifact. Build and Install: - stage: build + stage: Build script: - echo "**${CI_PROJECT_DIR}**" - git clone https://github.com/rest-for-physics/framework.git framework @@ -50,7 +50,7 @@ Build and Install: expire_in: 1 day Build Standalone: - stage: build-standalone + stage: Build Standalone script: - echo "**${CI_PROJECT_DIR}**" - cd ${CI_PROJECT_DIR}/ @@ -64,13 +64,13 @@ Build Standalone: - rm -rf ${CI_PROJECT_DIR}/install-standalone Load REST Libraries: - stage: loadRESTLibs + stage: Load REST Libraries script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - restRoot -b -q -01_NLDBD: - stage: examples +01.NLDBD: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/01.NLDBD/ @@ -78,8 +78,8 @@ Load REST Libraries: - geant4-config --version - restRoot -b -q Validate.C'("Run00001_NLDBD_Test.root")' -03_Fluorescence: - stage: examples +03.Fluorescence: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/03.Fluorescence/ @@ -88,8 +88,8 @@ Load REST Libraries: - geant4-config --version - restRoot -b -q Validate.C'("Run00111_g4Analysis_Gamma_rest.root")' -04_MuonScan: - stage: examples +04.MuonScan: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/04.MuonScan/ @@ -101,37 +101,37 @@ Load REST Libraries: - restG4 CosmicMuonsFromCircle.rml - restRoot -b -q ValidateCircle.C'("Run00111_restG4_MuonsFromCircle_rest.root")' -05_PandaXIII: - stage: examples +05.PandaXIII: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/05.PandaXIII/ - restG4 Xe136bb0n.rml - restRoot -b -q ValidateG4.C'("Xe136bb0n_n2E06.root")' -06_IonRecoils: - stage: examples +06.IonRecoils: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/06.IonRecoils/ - restG4 recoils.rml - restRoot -b -q Validate.C'("Run00001_F20_Recoils.root")' -07_FullChainDecay: - stage: examples +07.FullChainDecay: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/07.FullChainDecay/ - restG4 fullChain.rml - restG4 singleDecay.rml - - restRoot -b -q Validate.C'("Run00001_U238_FullChainDecay.root", 14)' + - restRoot -b -q Validate.C'("Run00001_U238_FullChainDecay.root", 16)' - restRoot -b -q Validate.C'("Run00002_U238_SingleChainDecay.root", 1)' - export REST_ISOTOPE="Be7" - restG4 singleDecay.rml - restRoot -b -q Validate.C'("Run00002_Be7_SingleChainDecay.root", 1)' -08_Alphas: - stage: examples +08.Alphas: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/08.Alphas/ @@ -151,28 +151,27 @@ Load REST Libraries: - restManager --c plots.rml --f "data/*g4Ana*root" --batch - restRoot -b -q Validate.C -09_Pb210_Shielding: - stage: examples +09.Pb210_Shielding: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/09.Pb210_Shield/ - restG4 Pb210.rml - restRoot -b -q Validate.C'("Run00001_Pb210_Shielding.root")' -10_Geometries: - stage: examples +10.Geometries: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/10.Geometries/ - restG4 Assembly.rml - restRoot -b -q Validate.C'("Run00001_Assembly_Assembly.root")' -11_Xrays: - stage: examples +11.Xrays: + stage: Run Examples script: - . ${CI_PROJECT_DIR}/install/thisREST.sh - cd ${CI_PROJECT_DIR}/install/examples/restG4/11.Xrays/ - restG4 xrays.rml - - restRoot -b -q Validate.C'("Run00001_Assembly_Assembly.root")' - restManager --c analysis.rml --f RestG4_xrays_00001_ArgonISO.root - restRoot -b -q GetQE.C'("Run_g4Analysis_00001_xrays.root")' diff --git a/examples/01.NLDBD/NLDBD.rml b/examples/01.NLDBD/NLDBD.rml index 57f5997f..52923f8c 100644 --- a/examples/01.NLDBD/NLDBD.rml +++ b/examples/01.NLDBD/NLDBD.rml @@ -44,7 +44,7 @@ By default REST units are mm, keV and degrees - + ///three types of source definition supported: diff --git a/examples/01.NLDBD/Validate.C b/examples/01.NLDBD/Validate.C index 0dc1d143..e4df5cfd 100644 --- a/examples/01.NLDBD/Validate.C +++ b/examples/01.NLDBD/Validate.C @@ -1,93 +1,126 @@ -Int_t Validate(string fname) { - TRestRun* run = new TRestRun(fname); +#include - if (run->GetParentRunNumber() != 0) { - cout << "Parent run number value : " << run->GetParentRunNumber() << endl; +Int_t Validate(const char* filename) { + cout << "Starting validation for '" << filename << "'" << endl; + + TRestRun run(filename); + + if (run.GetParentRunNumber() != 0) { + cout << "Parent run number value : " << run.GetParentRunNumber() << endl; cout << "The parent run number from restG4 generated file should be 0" << endl; return 1; } - if (run->GetRunNumber() != 1) { - cout << "Run number value : " << run->GetRunNumber() << endl; + if (run.GetRunNumber() != 1) { + cout << "Run number value : " << run.GetRunNumber() << endl; cout << "The run number on the validation chain should be 1 by default!" << endl; return 2; } - if (run->GetRunType() != "restG4") { - cout << "Run type : " << run->GetRunType() << endl; + if (run.GetRunType() != "restG4") { + cout << "Run type : " << run.GetRunType() << endl; cout << "The run type of restG4 generated data should be 'restG4'!" << endl; return 3; } - if (run->GetRunTag() != "NLDBD") { - cout << "Run tag : " << run->GetRunTag() << endl; + if (run.GetRunTag() != "NLDBD") { + cout << "Run tag : " << run.GetRunTag() << endl; cout << "The run tag of the basic validation test should be 'NLDBD'!" << endl; return 4; } - if (run->GetEntries() != 10) { - cout << "Run entries : " << run->GetEntries() << endl; - cout << "The NLDBD simulation is launched from gas. It should always generate 10 events." << endl; + if (run.GetEntries() != 100) { + cout << "Run entries : " << run.GetEntries() << endl; + cout << "The number of stored events should match the reference value of 100" << endl; return 5; } cout << "Testing reading of Geant4 metadata class" << endl; - TRestGeant4Metadata* g4Md = (TRestGeant4Metadata*)run->GetMetadataClass("TRestGeant4Metadata"); - if (!g4Md) { + TRestGeant4Metadata* geant4Metadata = (TRestGeant4Metadata*)run.GetMetadataClass("TRestGeant4Metadata"); + if (!geant4Metadata) { cout << "Problem reading Geant4 metadata class!" << endl; return 6; } - g4Md->PrintMetadata(); + geant4Metadata->PrintMetadata(); - if (g4Md->GetNumberOfActiveVolumes() != 1) { - cout << "The number of registered volumes is not 1!" << endl; + const bool isReferenceGeant4Version = geant4Metadata->GetGeant4Version() == "10.4.3"; + + if (geant4Metadata->GetNumberOfActiveVolumes() != 1) { + cout << "The number of registered does not match the reference value of 1" << endl; return 7; } - TRestGeant4Event* ev = (TRestGeant4Event*)run->GetInputEvent(); - run->GetEntry(9); + TRestGeant4Event* event = run.GetInputEvent(); + + double nEvents = run.GetEntries(); + + double averageTotalEnergy = 0; + const double averageTotalEnergyRef = (!isReferenceGeant4Version) ? 2280.96 : 2221.55; + + double averageSensitiveEnergy = 0; + const double averageSensitiveEnergyRef = (!isReferenceGeant4Version) ? 2280.96 : 2221.55; + + double averageNumberOfHits = 0; + const double averageNumberOfHitsRef = (!isReferenceGeant4Version) ? 3071.17 : 342.37; - cout << "Total energy : " << ev->GetTotalDepositedEnergy() << endl; - Int_t en = (Int_t)(100 * ev->GetTotalDepositedEnergy()); - cout << "Energy integer : " << en << endl; - cout << "Sensitive volume energy : " << ev->GetSensitiveVolumeEnergy() << endl; - cout << "Number of hits : " << ev->GetNumberOfHits() << endl; - cout << "Number of tracks : " << ev->GetNumberOfTracks() << endl; + double averageNumberOfTracks = 0; + const double averageNumberOfTracksRef = (!isReferenceGeant4Version) ? 2161.43 : 10.11; - Int_t X = (Int_t)(100 * ev->GetMeanPositionInVolume(0).X()); - Int_t Y = (Int_t)(100 * ev->GetMeanPositionInVolume(0).Y()); - Int_t Z = (Int_t)(100 * ev->GetMeanPositionInVolume(0).Z()); + TVector3 averagePosition = {}; + const TVector3 averagePositionRef = (!isReferenceGeant4Version) ? TVector3(-38.8987, 27.5536, 91.3969) + : TVector3(-17.8046, -32.5019, -31.8353); - cout << "x: " << X << " y: " << Y << " z: " << Z << endl; + const double tolerance = 0.001; - if (en != 245700) { - cout << "Error in total energy" << endl; + for (size_t i = 0; i < run.GetEntries(); i++) { + run.GetEntry(i); + + averageTotalEnergy += event->GetTotalDepositedEnergy() / nEvents; + averageSensitiveEnergy += event->GetSensitiveVolumeEnergy() / nEvents; + averageNumberOfHits += event->GetNumberOfHits() / nEvents; + averageNumberOfTracks += event->GetNumberOfTracks() / nEvents; + averagePosition += event->GetMeanPositionInVolume(0) * (1.0 / nEvents); + } + + cout << "Average total energy: " << averageTotalEnergy << " keV" << endl; + cout << "Average sensitive energy: " << averageSensitiveEnergy << " keV" << endl; + cout << "Average number of hits: " << averageNumberOfHits << endl; + cout << "Average number of tracks: " << averageNumberOfTracks << endl; + cout << "Average position: (" << averagePosition.x() << ", " << averagePosition.y() << ", " + << averagePosition.z() << ") mm" << endl; + + if (TMath::Abs(averageNumberOfHits - averageNumberOfHitsRef) / averageNumberOfHitsRef > tolerance) { + cout << "The average number of hits does not match the reference value of " << averageNumberOfHitsRef + << endl; return 8; } - //en = (Int_t)(100 * ev->GetSensitiveVolumeEnergy()); - //if (en != 245699) { - // cout << "Error in total energy" << endl; - // return 9; - //} + if (TMath::Abs(averageNumberOfTracks - averageNumberOfTracksRef) / averageNumberOfTracksRef > tolerance) { + cout << "The average number of tracks does not match the reference value of " + << averageNumberOfTracksRef << endl; + return 9; + } - if (ev->GetNumberOfHits() != 401) { - cout << "Error in the number of hits" << endl; + if (TMath::Abs(averageSensitiveEnergy - averageSensitiveEnergyRef) / averageSensitiveEnergyRef > + tolerance) { + cout << "The average sensitive volume energy does not match the reference value of " + << averageSensitiveEnergyRef << endl; return 10; } - if (ev->GetNumberOfTracks() != 11) { - cout << "Error in the number of tracks" << endl; + if (TMath::Abs(averageTotalEnergy - averageTotalEnergyRef) / averageTotalEnergyRef > tolerance) { + cout << "The average total energy does not match the reference value of " << averageTotalEnergyRef + << endl; return 11; } - if (X != -52709 || Y != 37592 || Z != 59585) { - cout << "Error in the event mean position" << endl; + + if (TMath::Abs(averagePosition.Mag() - averagePositionRef.Mag()) / averagePositionRef.Mag() > tolerance) { + cout << "The average position does not match the reference value of " + << "(" << averagePositionRef.x() << ", " << averagePositionRef.y() << ", " + << averagePositionRef.z() << ") mm" << endl; return 12; } - cout << "All tests passed! [\033[32mOK\033[0m]\n"; - // Other tests like opening other metadata classes. Detector TGeoManager, etc. - return 0; } diff --git a/examples/04.MuonScan/CosmicMuonsFromCircle.rml b/examples/04.MuonScan/CosmicMuonsFromCircle.rml index c6ab44d1..59e0a853 100644 --- a/examples/04.MuonScan/CosmicMuonsFromCircle.rml +++ b/examples/04.MuonScan/CosmicMuonsFromCircle.rml @@ -65,7 +65,6 @@ - diff --git a/examples/04.MuonScan/CosmicMuonsFromWall.rml b/examples/04.MuonScan/CosmicMuonsFromWall.rml index 46ed1c42..a88d5574 100644 --- a/examples/04.MuonScan/CosmicMuonsFromWall.rml +++ b/examples/04.MuonScan/CosmicMuonsFromWall.rml @@ -52,24 +52,23 @@ First concept author G. Luzón (16-11-2020) --> - + - - + + - - - - - + + + + diff --git a/examples/04.MuonScan/Muon.rml b/examples/04.MuonScan/Muon.rml index f13c0b6e..da6315d7 100644 --- a/examples/04.MuonScan/Muon.rml +++ b/examples/04.MuonScan/Muon.rml @@ -17,7 +17,7 @@ - + @@ -48,24 +48,23 @@ - + - - + + - - - - - + + + + diff --git a/examples/04.MuonScan/Validate.C b/examples/04.MuonScan/Validate.C index 205f12a2..c4ea54f0 100644 --- a/examples/04.MuonScan/Validate.C +++ b/examples/04.MuonScan/Validate.C @@ -1,103 +1,140 @@ -Int_t Validate(string fname) { - gSystem->Load("libRestFramework.so"); - gSystem->Load("libRestGeant4.so"); +#include - TRestRun* run = new TRestRun(fname); +Int_t Validate(const char* filename) { + cout << "Starting validation for '" << filename << "'" << endl; - if (run->GetParentRunNumber() != 0) { - cout << "Parent run number value : " << run->GetParentRunNumber() << endl; + TRestRun run(filename); + + if (run.GetParentRunNumber() != 0) { + cout << "Parent run number value : " << run.GetParentRunNumber() << endl; cout << "The parent run number from restG4 generated file should be 0" << endl; return 1; } - if (run->GetRunNumber() != 111) { - cout << "Run number value : " << run->GetRunNumber() << endl; + if (run.GetRunNumber() != 111) { + cout << "Run number value : " << run.GetRunNumber() << endl; cout << "The run number on the Muon example should be 111 by default!" << endl; return 2; } - if (run->GetRunType() != "restG4") { - cout << "Run type : " << run->GetRunType() << endl; + if (run.GetRunType() != "restG4") { + cout << "Run type : " << run.GetRunType() << endl; cout << "The run type of restG4 generated data should be 'restG4'!" << endl; return 3; } - if (run->GetRunTag() != "MuonsFromPoint") { - cout << "Run tag : " << run->GetRunTag() << endl; + if (run.GetRunTag() != "MuonsFromPoint") { + cout << "Run tag : " << run.GetRunTag() << endl; cout << "The run tag of the basic validation test should be 'MuonsFromPoint'!" << endl; return 4; } - if (run->GetEntries() != 10) { - cout << "Run entries : " << run->GetEntries() << endl; - cout << "The Muon simulation should always generate 10 events." << endl; + if (run.GetEntries() != 1000) { + cout << "Run entries : " << run.GetEntries() << endl; + cout << "The number of stored events should match the reference value of 1000" << endl; return 5; } cout << "Testing reading of Geant4 metadata class" << endl; - TRestGeant4Metadata* g4Md = (TRestGeant4Metadata*)run->GetMetadataClass("TRestGeant4Metadata"); - if (!g4Md) { + TRestGeant4Metadata* geant4Metadata = (TRestGeant4Metadata*)run.GetMetadataClass("TRestGeant4Metadata"); + if (!geant4Metadata) { cout << "Problem reading Geant4 metadata class!" << endl; return 6; } - g4Md->PrintMetadata(); + geant4Metadata->PrintMetadata(); - if (g4Md->GetNumberOfActiveVolumes() != 2) { - cout << "The number of registered volumes is not 2!" << endl; + if (geant4Metadata->GetNumberOfActiveVolumes() != 2) { + cout << "The number of registered does not match the reference value of 2" << endl; return 7; } - TRestGeant4Event* ev = (TRestGeant4Event*)run->GetInputEvent(); - run->GetEntry(9); + TRestGeant4Event* event = run.GetInputEvent(); + + double nEvents = run.GetEntries(); + + double averageTotalEnergy = 0; + constexpr double averageTotalEnergyRef = 18.0321; + + double averageSensitiveEnergy = 0; + constexpr double averageSensitiveEnergyRef = 8.78066; + + double averageNumberOfTracks = 0; + constexpr double averageNumberOfTracksRef = 4.924; + + double averageNumberOfHitsVolume0 = 0; + constexpr double averageNumberOfHitsVolume0Ref = 60.712; - cout << "Total energy : " << ev->GetTotalDepositedEnergy() << endl; - Int_t en = (Int_t)(100 * ev->GetTotalDepositedEnergy()); - cout << "Energy integer : " << en << endl; - cout << "Sensitive volume energy : " << ev->GetSensitiveVolumeEnergy() << endl; - cout << "Number of hits : " << ev->GetNumberOfHits() << endl; - cout << "Number of tracks : " << ev->GetNumberOfTracks() << endl; + double averageNumberOfHitsVolume1 = 0; + constexpr double averageNumberOfHitsVolume1Ref = 68.061; - Int_t X = (Int_t)(100 * ev->GetMeanPositionInVolume(0).X()); - Int_t Y = (Int_t)(100 * ev->GetMeanPositionInVolume(0).Y()); - Int_t Z = (Int_t)(100 * ev->GetMeanPositionInVolume(0).Z()); + TVector3 averagePosition = {}; + const TVector3 averagePositionRef = {-0.0305454, -0.384008, -299.54}; - cout << "x: " << X << " y: " << Y << " z: " << Z << endl; + constexpr double tolerance = 0.001; - if (en != 1206) { - cout << "Error in total energy" << endl; + for (size_t i = 0; i < run.GetEntries(); i++) { + run.GetEntry(i); + + averageTotalEnergy += event->GetTotalDepositedEnergy() / nEvents; + averageSensitiveEnergy += event->GetSensitiveVolumeEnergy() / nEvents; + averageNumberOfHitsVolume0 += event->GetNumberOfHits(0) / nEvents; + averageNumberOfHitsVolume1 += event->GetNumberOfHits(1) / nEvents; + averageNumberOfTracks += event->GetNumberOfTracks() / nEvents; + averagePosition += event->GetMeanPositionInVolume(0) * (1.0 / nEvents); + } + + cout << "Average total energy: " << averageTotalEnergy << " keV" << endl; + cout << "Average sensitive energy: " << averageSensitiveEnergy << " keV" << endl; + cout << "Average number of hits in volume 0: " << averageNumberOfHitsVolume0 << endl; + cout << "Average number of hits in volume 1: " << averageNumberOfHitsVolume1 << endl; + cout << "Average number of tracks: " << averageNumberOfTracks << endl; + cout << "Average position: (" << averagePosition.x() << ", " << averagePosition.y() << ", " + << averagePosition.z() << ") mm" << endl; + + if (TMath::Abs(averageNumberOfTracks - averageNumberOfTracksRef) / averageNumberOfTracksRef > + tolerance) { + cout << "The average number of tracks does not match the reference value of " + << averageNumberOfTracksRef << endl; return 8; } - en = (Int_t)(100 * ev->GetSensitiveVolumeEnergy()); - if (en != 496) { - cout << "Error in total energy" << endl; + if (TMath::Abs(averageNumberOfHitsVolume0 - averageNumberOfHitsVolume0Ref) / + averageNumberOfHitsVolume0Ref > + tolerance) { + cout << "The average number of hits in volume 0 does not match the reference value of " + << averageNumberOfHitsVolume0Ref << endl; return 9; } - if (ev->GetNumberOfHits() != 102) { - cout << "Error in the number of hits" << endl; + if (TMath::Abs(averageNumberOfHitsVolume1 - averageNumberOfHitsVolume1Ref) / + averageNumberOfHitsVolume1Ref > + tolerance) { + cout << "The average number of hits in volume 1 does not match the reference value of " + << averageNumberOfHitsVolume1Ref << endl; return 10; } - if (ev->GetNumberOfTracks() != 1) { - cout << "Error in the number of tracks" << endl; + if (TMath::Abs(averageSensitiveEnergy - averageSensitiveEnergyRef) / averageSensitiveEnergyRef > + tolerance) { + cout << "The average sensitive volume energy does not match the reference value of " + << averageSensitiveEnergyRef << endl; return 11; } - if (X != 0 || Y != 0 || Z != -30180) { - cout << "Error in the event mean position" << endl; + + if (TMath::Abs(averageTotalEnergy - averageTotalEnergyRef) / averageTotalEnergyRef > tolerance) { + cout << "The average total energy does not match the reference value of " << averageTotalEnergyRef + << endl; return 12; } - cout << "Number of hits in active volume 0: " << ev->GetNumberOfHits(0) << endl; - cout << "Number of hits in active volume 1: " << ev->GetNumberOfHits(1) << endl; - if (ev->GetNumberOfHits(0) != 51 || ev->GetNumberOfHits(1) != 51) { - cout << "Error in number of hits at the active volumes" << endl; + if (TMath::Abs(averagePosition.Mag() - averagePositionRef.Mag()) / averagePositionRef.Mag() > + tolerance) { + cout << "The average position does not match the reference value of " + << "(" << averagePositionRef.x() << ", " << averagePositionRef.y() << ", " + << averagePositionRef.z() << ") mm" << endl; return 13; } - cout << "All tests passed! [\033[32mOK\033[0m]\n"; - // Other tests like opening other metadata classes. Detector TGeoManager, etc. - return 0; } diff --git a/examples/04.MuonScan/ValidateCircle.C b/examples/04.MuonScan/ValidateCircle.C index 63d9b4e0..2a353bfc 100644 --- a/examples/04.MuonScan/ValidateCircle.C +++ b/examples/04.MuonScan/ValidateCircle.C @@ -1,15 +1,15 @@ +#include -Int_t ValidateCircle(string fname) { - gSystem->Load("/usr/local/rest-for-physics/lib/libRestFramework.so"); - gSystem->Load("/usr/local/rest-for-physics/lib/libRestGeant4.so"); +Int_t ValidateCircle(const char* filename) { + cout << "Starting validation for '" << filename << "'" << endl; - TRestRun run(fname); + TRestRun run(filename); TRestGeant4Event* event = run.GetInputEvent(); if (run.GetRunTag() != "MuonsFromCircle") { cout << "Run tag : " << run.GetRunTag() << endl; cout << "The run tag of the basic validation test should be 'MuonsFromCircle'" << endl; - return 4; + return 1; } Double_t rMean = 0; @@ -34,26 +34,26 @@ Int_t ValidateCircle(string fname) { if (rMean < 2.75 || rMean > 3.25) { cout << "The average radius of the distribution is wrong!" << endl; cout << "R_mean (cm): " << rMean << endl; - return 5; + return 2; } if (rMin > 0.5) { cout << "The minimum radius of the distribution is wrong!" << endl; cout << "R_min (cm): " << rMin << endl; - return 6; + return 3; } if (rMax > 4.0 || rMax < 3.75) { cout << "The maximum radius of the distribution is wrong!" << endl; cout << "R_max (cm): " << rMax << endl; - return 7; + return 4; } cout << "Run entries: " << run.GetEntries() << endl; if (run.GetEntries() < 600 || run.GetEntries() > 750) { cout << "The number of entries is wrong!" << endl; cout << "Number of entries : " << run.GetEntries() << endl; - return 8; + return 5; } cout << "All tests passed! [\033[32mOK\033[0m]\n"; diff --git a/examples/04.MuonScan/ValidateWall.C b/examples/04.MuonScan/ValidateWall.C index e9da9f15..7e6da002 100644 --- a/examples/04.MuonScan/ValidateWall.C +++ b/examples/04.MuonScan/ValidateWall.C @@ -1,38 +1,38 @@ +#include -Int_t ValidateWall(string fname) { - gSystem->Load("libRestFramework.so"); - gSystem->Load("libRestGeant4.so"); +Int_t ValidateWall(const char* filename) { + cout << "Starting validation for '" << filename << "'" << endl; - TRestRun* run = new TRestRun(fname); - TRestGeant4Event* ev = (TRestGeant4Event*)run->GetInputEvent(); + TRestRun run(filename); + TRestGeant4Event* event = run.GetInputEvent(); - if (run->GetRunTag() != "MuonsFromWall") { - cout << "Run tag : " << run->GetRunTag() << endl; + if (run.GetRunTag() != "MuonsFromWall") { + cout << "Run tag : " << run.GetRunTag() << endl; cout << "The run tag of the basic validation test should be 'MuonsFromWall!" << endl; - return 4; + return 1; } Double_t r = 0; - for (Int_t n = 0; n < run->GetEntries(); n++) { - run->GetEntry(n); - Double_t x = ev->GetPrimaryEventOrigin().X(); - Double_t y = ev->GetPrimaryEventOrigin().Y(); + for (Int_t n = 0; n < run.GetEntries(); n++) { + run.GetEntry(n); + Double_t x = event->GetPrimaryEventOrigin().X(); + Double_t y = event->GetPrimaryEventOrigin().Y(); r += x * x + y * y; } - r /= run->GetEntries(); + r /= run.GetEntries(); if ((Int_t)(1000. * r) / 1000 < 10000 || (Int_t)(1000. * r) / 1000 > 20000) { cout << "The average radius of the distribution is wrong!" << endl; cout << "R: " << (Int_t)(1000. * r) / 1000 << endl; - return 5; + return 2; } - cout << "Run entries: " << run->GetEntries() << endl; - if (run->GetEntries() < 350 || run->GetEntries() > 450) { + cout << "Run entries: " << run.GetEntries() << endl; + if (run.GetEntries() < 350 || run.GetEntries() > 450) { cout << "The number of entries is not between 350 and 450!" << endl; - cout << "Number of entries : " << run->GetEntries() << endl; - return 6; + cout << "Number of entries : " << run.GetEntries() << endl; + return 3; } cout << "All tests passed! [\033[32mOK\033[0m]\n"; diff --git a/examples/07.FullChainDecay/Validate.C b/examples/07.FullChainDecay/Validate.C index 7daac3b3..d01296ec 100644 --- a/examples/07.FullChainDecay/Validate.C +++ b/examples/07.FullChainDecay/Validate.C @@ -1,31 +1,28 @@ -Int_t Validate(string fname, Int_t nDaughters) { +Int_t Validate(const char* filename, Int_t nDaughters) { gSystem->Load("libRestFramework.so"); gSystem->Load("libRestGeant4.so"); - TRestRun* run = new TRestRun(fname); + TRestRun run(filename); - TRestGeant4Event* ev = (TRestGeant4Event*)run->GetInputEvent(); + TRestGeant4Event* event = run.GetInputEvent(); - std::vector evTag; - for (int n = 0; n < run->GetEntries(); n++) { - run->GetEntry(n); - evTag.push_back((string)ev->GetSubEventTag()); + std::set eventTagsUnique; + for (int n = 0; n < run.GetEntries(); n++) { + run.GetEntry(n); + eventTagsUnique.insert((string)event->GetSubEventTag()); } - std::sort(evTag.begin(), evTag.end()); - - auto iter = std::unique(evTag.begin(), evTag.end()); - - evTag.erase(iter, evTag.end()); - - cout << "Daughter isotopes: " << evTag.size() << endl; - for (unsigned int n = 0; n < evTag.size(); n++) cout << evTag[n] << " "; + cout << "Daughter isotopes: " << eventTagsUnique.size() << endl; + for (const auto& tag : eventTagsUnique) { + cout << tag << " "; + } cout << endl; - if (evTag.size() != nDaughters) { - cout << "Wrong number of isotopes found!" << endl; - return 13; + if (eventTagsUnique.size() != nDaughters) { + cout << "Wrong number of isotopes found! " << eventTagsUnique.size() << " vs. reference value of " + << nDaughters << endl; + return 1; } cout << "All tests passed! [\033[32mOK\033[0m]\n"; diff --git a/examples/07.FullChainDecay/fullChain.rml b/examples/07.FullChainDecay/fullChain.rml index 9b5c98eb..47b36846 100644 --- a/examples/07.FullChainDecay/fullChain.rml +++ b/examples/07.FullChainDecay/fullChain.rml @@ -37,7 +37,7 @@ By default REST units are mm, keV and degrees - + diff --git a/examples/07.FullChainDecay/singleDecay.rml b/examples/07.FullChainDecay/singleDecay.rml index 2dea07eb..c91cc5dd 100644 --- a/examples/07.FullChainDecay/singleDecay.rml +++ b/examples/07.FullChainDecay/singleDecay.rml @@ -37,7 +37,7 @@ By default REST units are mm, keV and degrees - + diff --git a/examples/08.Alphas/Validate.C b/examples/08.Alphas/Validate.C index be20bc98..caf5fed9 100644 --- a/examples/08.Alphas/Validate.C +++ b/examples/08.Alphas/Validate.C @@ -1,75 +1,76 @@ Int_t Validate() { - TFile* f = new TFile("plots.root"); + TFile* file = TFile::Open("plots.root"); - TH1F* h1 = (TH1F*)f->Get("1MeV_5um"); - TH1F* h2 = (TH1F*)f->Get("5MeV_5um"); - TH1F* h3 = (TH1F*)f->Get("5MeV_1um"); + TH1F* h1 = file->Get("1MeV_5um"); + TH1F* h2 = file->Get("5MeV_5um"); + TH1F* h3 = file->Get("5MeV_1um"); cout << "Entries: h1: " << h1->Integral() << ", h2: " << h2->Integral() << ", h3: " << h3->Integral() << endl; + TRestRun run("data/Run_g4Analysis_1MeV_5um.root"); + TRestAnalysisTree* analysisTree = run.GetAnalysisTree(); - TRestRun* run = new TRestRun("data/Run_g4Analysis_1MeV_5um.root"); - TRestAnalysisTree* aT = run->GetAnalysisTree(); + run.GetEntry(100); - run->GetEntry(100); + Double_t thetaSample = analysisTree->GetObservableValue("g4Ana_thetaPrimary"); + constexpr Double_t thetaSampleRef = 1.77307; - Int_t thetaInt = (int)(1000 * aT->GetObservableValue("g4Ana_thetaPrimary")); - Int_t phiInt = (int)(1000 * aT->GetObservableValue("g4Ana_phiPrimary")); - Double_t theta = aT->GetObservableAverage("g4Ana_thetaPrimary"); - Double_t phi = aT->GetObservableAverage("g4Ana_phiPrimary"); + Double_t phiSample = analysisTree->GetObservableValue("g4Ana_phiPrimary"); + constexpr Double_t phiSampleRef = -1.59582; - cout << "entry 100, thetaInt: " << thetaInt << ", phiInt: " << phiInt << endl; - cout << "average theta: " << theta << ", average phi: " << phi << endl; + Double_t thetaAverage = analysisTree->GetObservableAverage("g4Ana_thetaPrimary"); + constexpr Double_t thetaAverageRef = 1.5767; - - if (h1->Integral() != 1569) { + Double_t phiAverage = analysisTree->GetObservableAverage("g4Ana_phiPrimary"); + constexpr Double_t phiAverageRef = 0.0727066; + + cout << "entry 100, theta: " << thetaSample << ", phi: " << phiSample << endl; + cout << "average theta: " << thetaAverage << ", average phi: " << phiAverage << endl; + + if (h1->Integral() != 1960) { cout << "Wrong number of alphas produced at 1MeV_5um!" << endl; cout << "Histogram contains : " << h1->Integral() << endl; return 10; } - if (h2->Integral() != 7286) { + if (h2->Integral() != 7623) { cout << "Wrong number of alphas produced at 5MeV_5um!" << endl; cout << "Histogram contains : " << h2->Integral() << endl; return 20; } - if (h3->Integral() != 9401) { + if (h3->Integral() != 9554) { cout << "Wrong number of alphas produced at 5MeV_1um!" << endl; cout << "Histogram contains : " << h3->Integral() << endl; return 30; } - if (thetaInt != 2787) { + if (TMath::Abs(thetaSample - thetaSampleRef) / thetaSampleRef >= 0.01) { cout << "Wrong theta angle value for entry 100!" << endl; - cout << "Theta value is : " << thetaInt << endl; + cout << "Theta value is : " << thetaSample << endl; return 60; } - - if (phiInt != 1855) { + if (TMath::Abs(phiSample - phiSampleRef) / phiSampleRef >= 0.01) { cout << "Wrong phi angle value for entry 100!" << endl; - cout << "Phi value is : " << phiInt << endl; + cout << "Phi value is : " << phiSample << endl; return 60; } - - if ((int)(1000 * theta) != 1566) { + if (TMath::Abs(thetaAverage - thetaAverageRef) / thetaAverageRef >= 0.01) { cout << "Wrong theta angle average!" << endl; - cout << "Theta angle average : " << theta << " while it should be : " << 1.589 << endl; + cout << "Theta angle average : " << thetaAverage << " while it should be : " << 1.589 << endl; return 80; } - - if ((int)(1000 * phi) != -4) { + if (TMath::Abs(phiAverage - phiAverageRef) / phiAverageRef >= 0.01) { cout << "Wrong phi angle average!" << endl; - cout << "Phi angle average : " << phi << " while it should be : " << 0.027 << endl; + cout << "Phi angle average : " << phiAverage << " while it should be : " << 0.027 << endl; return 90; } - cout << "All tests passed! [\033[32mOK\033[0m]\n"; return 0; } diff --git a/examples/08.Alphas/alphas.rml b/examples/08.Alphas/alphas.rml index 47c03ab9..c70cbf37 100644 --- a/examples/08.Alphas/alphas.rml +++ b/examples/08.Alphas/alphas.rml @@ -20,15 +20,15 @@ By default REST units are mm, keV and degrees - - + + - + @@ -40,15 +40,15 @@ By default REST units are mm, keV and degrees - + - + - + @@ -66,24 +66,24 @@ By default REST units are mm, keV and degrees // Use only one EM physics list - + - + - - + + - - - - - + + + + + diff --git a/examples/08.Alphas/analysis.rml b/examples/08.Alphas/analysis.rml index 1c581e0a..7ccfa48e 100644 --- a/examples/08.Alphas/analysis.rml +++ b/examples/08.Alphas/analysis.rml @@ -24,7 +24,7 @@ internal values are saved. - + @@ -37,12 +37,12 @@ internal values are saved. // observable = all will add all NON `custom` observables - - - - - - + + + + + + diff --git a/examples/08.Alphas/plots.rml b/examples/08.Alphas/plots.rml index 96abfff9..402101ad 100644 --- a/examples/08.Alphas/plots.rml +++ b/examples/08.Alphas/plots.rml @@ -1,37 +1,35 @@ - + - + - + - + - + - - - - + + + + - - - - + + + + - - - - + + + + - + - + - + - - - diff --git a/examples/09.Pb210_Shield/Pb210.rml b/examples/09.Pb210_Shield/Pb210.rml index d25e5f7e..0f06cb8e 100644 --- a/examples/09.Pb210_Shield/Pb210.rml +++ b/examples/09.Pb210_Shield/Pb210.rml @@ -36,7 +36,7 @@ By default REST units are mm, keV and degrees - + diff --git a/examples/09.Pb210_Shield/Validate.C b/examples/09.Pb210_Shield/Validate.C index 00f6c830..2515c3f3 100644 --- a/examples/09.Pb210_Shield/Validate.C +++ b/examples/09.Pb210_Shield/Validate.C @@ -1,14 +1,17 @@ -Int_t Validate(string fname) { +Int_t Validate(const char* filename) { gSystem->Load("libRestFramework.so"); gSystem->Load("libRestGeant4.so"); - TRestRun* run = new TRestRun(fname); + TRestRun run(filename); - if (run->GetEntries() != 28) { - cout << "Entries: " << run->GetEntries() << endl; - cout << "There was a problem simulation Pb210. The number of entries should be 28" << endl; - return 10; + constexpr double tolerance = 0.001; + constexpr int numberOfEntriesRef = 252; + + if (TMath::Abs(run.GetEntries() - numberOfEntriesRef) / numberOfEntriesRef > tolerance) { + cout << "Number of entries: " << run.GetEntries() << "does not match the reference value of " + << numberOfEntriesRef << endl; + return 1; } cout << "All tests passed! [\033[32mOK\033[0m]\n"; diff --git a/examples/11.Xrays/GetQE.C b/examples/11.Xrays/GetQE.C index f8f6ae36..070e5fd0 100644 --- a/examples/11.Xrays/GetQE.C +++ b/examples/11.Xrays/GetQE.C @@ -1,85 +1,85 @@ #include "TRestRun.h" -Int_t GetQE(string inputFile) { - TRestRun* run = new TRestRun(inputFile); - TRestAnalysisTree* aTree = run->GetAnalysisTree(); +Int_t GetQE(const char* inputFile) { + TRestRun run(inputFile); + TRestAnalysisTree* aTree = run.GetAnalysisTree(); - TH1D* edepHist = new TH1D("edepHist", "Energy deposited in detector", 150, 0, 15); - TH1D* primaryHist = new TH1D("primaryHist", "Energy of generated primary", 150, 0, 15); + auto eDepHist = new TH1D("eDepHist", "Energy deposited in detector", 150, 0, 15); + auto primaryHist = new TH1D("primaryHist", "Energy of generated primary", 150, 0, 15); - for (int n = 0; n < run->GetEntries(); n++) { - run->GetEntry(n); + for (int n = 0; n < run.GetEntries(); n++) { + run.GetEntry(n); Double_t eDep = aTree->GetObservableValue("g4Ana_totalEdep"); Double_t ePrimary = aTree->GetObservableValue("g4Ana_energyPrimary"); - if (eDep > 0) edepHist->Fill(eDep); + if (eDep > 0) eDepHist->Fill(eDep); primaryHist->Fill(ePrimary); } FILE* f = fopen("output.txt", "wt"); - for (int n = 0; n < edepHist->GetNbinsX(); n++) { - Double_t edepCounts = edepHist->GetBinContent(n + 1); + for (int n = 0; n < eDepHist->GetNbinsX(); n++) { + Double_t eDepCounts = eDepHist->GetBinContent(n + 1); Double_t primaryCounts = primaryHist->GetBinContent(n + 1); - Double_t eff = edepCounts / primaryCounts; - Double_t err = TMath::Sqrt(edepCounts) / primaryCounts; - fprintf(f, "%4.2f\t%5.0lf\t%5.0lf\t%4.4lf\t%4.4lf\n", edepHist->GetBinCenter(n + 1), edepCounts, + Double_t eff = eDepCounts / primaryCounts; + Double_t err = TMath::Sqrt(eDepCounts) / primaryCounts; + fprintf(f, "%4.2f\t%5.0lf\t%5.0lf\t%4.4lf\t%4.4lf\n", eDepHist->GetBinCenter(n + 1), eDepCounts, primaryCounts, eff, err); } fclose(f); // Plot - TCanvas* c1 = new TCanvas("c1", "response matrix", 1200, 800); + auto c1 = new TCanvas("c1", "response matrix", 1200, 800); c1->GetFrame()->SetBorderSize(6); c1->GetFrame()->SetBorderMode(-1); - TPad* pad1 = new TPad("pad1", "This is pad1", 0.01, 0.02, 0.99, 0.97); - pad1->Divide(2, 2); - pad1->Draw(); + auto pad = new TPad("pad", "This is pad", 0.01, 0.02, 0.99, 0.97); + pad->Divide(2, 2); + pad->Draw(); - TH2D* matrix_hist = new TH2D("matrix hist", "response matrix", 150, 0, 15, 150, 0, 15); - matrix_hist->GetXaxis()->SetTitle("deposited Energy [keV]"); - matrix_hist->GetYaxis()->SetTitle("primary Energy [keV]"); + auto matrixHist = new TH2D("matrix hist", "response matrix", 150, 0, 15, 150, 0, 15); + matrixHist->GetXaxis()->SetTitle("deposited Energy [keV]"); + matrixHist->GetYaxis()->SetTitle("primary Energy [keV]"); gPad->SetLogz(); // gPad->SetTheta(90); // gPad->SetPhi(0); // gPad->SetLeftMargin(0.15); // gPad->SetBottomMargin(0.15); - for (int n = 0; n < run->GetEntries(); n++) { - run->GetEntry(n); + for (int n = 0; n < run.GetEntries(); n++) { + run.GetEntry(n); Double_t eDep = aTree->GetObservableValue("g4Ana_totalEdep"); Double_t ePrimary = aTree->GetObservableValue("g4Ana_energyPrimary"); - if (eDep > 0) matrix_hist->Fill(eDep, ePrimary); + if (eDep > 0) matrixHist->Fill(eDep, ePrimary); } - TH1D* prHist = matrix_hist->ProjectionY("PrimaryEnergy"); - pad1->cd(1); + TH1D* prHist = matrixHist->ProjectionY("PrimaryEnergy"); + pad->cd(1); prHist->Draw(); - TH1D* depHist = matrix_hist->ProjectionX("DepositEnergy"); - pad1->cd(2); + TH1D* depHist = matrixHist->ProjectionX("DepositEnergy"); + pad->cd(2); depHist->Draw(); - pad1->cd(3); + pad->cd(3); primaryHist->Draw(); - pad1->cd(4); - edepHist->Draw(); + pad->cd(4); + eDepHist->Draw(); - // gPad->cd(3); - // matrix_hist->Draw("SURF1"); + // gPad->cd(3); + // matrixHist->Draw("SURF1"); c1->Print("montecarlo.png"); - TH1D* effHisto = new TH1D("Efficiency", "Argon isobutane at 1 bar", 150, 0, 15); + auto effHisto = new TH1D("Efficiency", "Argon isobutane at 1 bar", 150, 0, 15); for (int n = 1; n <= effHisto->GetNbinsX(); n++) - effHisto->SetBinContent(n, edepHist->GetBinContent(n) / primaryHist->GetBinContent(n)); + effHisto->SetBinContent(n, eDepHist->GetBinContent(n) / primaryHist->GetBinContent(n)); - TCanvas* c2 = new TCanvas("c2", "Efficiency", 800, 600); + auto c2 = new TCanvas("c2", "Efficiency", 800, 600); c2->GetFrame()->SetBorderSize(6); c2->GetFrame()->SetBorderMode(-1); @@ -88,9 +88,11 @@ Int_t GetQE(string inputFile) { c2->Print("efficiency.png"); - cout << "Integral: " << edepHist->Integral() << endl; - if (edepHist->Integral() < 500) { - cout << "Something went wrong. Number of counts inside deposits integral too low" << endl; + constexpr double eDepIntegralReference = 500.0; + cout << "Integral: " << eDepHist->Integral() << endl; + if (eDepHist->Integral() < eDepIntegralReference) { + cout << "Number of counts inside deposits integral " << eDepHist->Integral() + << " does not match reference value of " << eDepIntegralReference << endl; return 1; } diff --git a/examples/11.Xrays/README.md b/examples/11.Xrays/README.md index b3566ec6..80c4f95b 100644 --- a/examples/11.Xrays/README.md +++ b/examples/11.Xrays/README.md @@ -2,7 +2,7 @@ **Author:** Javier Galan (javier.galan@unizar.es) -This example uses a x-ray gun with photons in the (0,15)keV range hitting a gas volume filled with Xenon in order to generate an efficiency response matrix. +This example uses an x-ray gun with photons in the (0,15)keV range hitting a gas volume filled with Xenon in order to generate an efficiency response matrix. #### Detector geometry diff --git a/examples/11.Xrays/analysis.rml b/examples/11.Xrays/analysis.rml index cd6d6a7f..c4074e69 100644 --- a/examples/11.Xrays/analysis.rml +++ b/examples/11.Xrays/analysis.rml @@ -24,7 +24,7 @@ internal values are saved. - + @@ -37,8 +37,8 @@ internal values are saved. // observable = all will add all NON `custom` observables - - + + diff --git a/examples/11.Xrays/xrays.rml b/examples/11.Xrays/xrays.rml index 859ecdfd..ea10513e 100644 --- a/examples/11.Xrays/xrays.rml +++ b/examples/11.Xrays/xrays.rml @@ -67,23 +67,22 @@ By default REST units are mm, keV and degrees - + - - + + - - - - - + + + + + - diff --git a/include/EventAction.h b/include/EventAction.h index 4b36c6e1..c787661e 100644 --- a/include/EventAction.h +++ b/include/EventAction.h @@ -6,6 +6,7 @@ #include #include #include +#include #include #include @@ -20,6 +21,8 @@ class EventAction : public G4UserEventAction { virtual void EndOfEventAction(const G4Event*); private: + TStopwatch fTimer; + Double_t absDouble(Double_t x) { if (x > 0) return x; return -x; diff --git a/restG4.cxx b/restG4.cxx index 254e6b5a..e8c6c29d 100644 --- a/restG4.cxx +++ b/restG4.cxx @@ -60,7 +60,7 @@ TH2D* spatialDistribution[maxBiasingVolumes]; TH1D initialEnergySpectrum; TH1D initialAngularDistribution; -Int_t N_events; +Int_t nEvents; int main(int argc, char** argv) { auto start_time = chrono::steady_clock::now(); @@ -74,6 +74,12 @@ int main(int argc, char** argv) { /// Separating relative path and pure RML filename char* inputConfigFile = const_cast(commandLineParameters.rmlFile.Data()); + + if (!TRestTools::CheckFileIsAccessible(inputConfigFile)) { + cout << "Input rml file: " << inputConfigFile << " not found, please check file name" << endl; + exit(1); + } + const auto [inputRmlPath, inputRmlClean] = TRestTools::SeparatePathAndName(inputConfigFile); if (!filesystem::path(inputRmlPath).empty()) { @@ -81,14 +87,12 @@ int main(int argc, char** argv) { } restG4Metadata = new TRestGeant4Metadata(inputRmlClean.c_str()); + restG4Metadata->SetGeant4Version(TRestTools::Execute("geant4-config --version")); if (!commandLineParameters.geometryFile.IsNull()) { restG4Metadata->SetGdmlFilename(commandLineParameters.geometryFile.Data()); } - string geant4Version = TRestTools::Execute("geant4-config --version"); - restG4Metadata->SetGeant4Version(geant4Version); - // We need to process and generate a new GDML for several reasons. // 1. ROOT6 has problem loading math expressions in gdml file // 2. We allow file entities to be http remote files @@ -250,7 +254,7 @@ int main(int argc, char** argv) { visManager->Initialize(); #endif - N_events = restG4Metadata->GetNumberOfEvents(); + nEvents = restG4Metadata->GetNumberOfEvents(); // We pass the volume definition to Stepping action so that it records gammas // entering in We pass also the biasing spectrum so that gammas energies // entering the volume are recorded @@ -264,8 +268,8 @@ int main(int argc, char** argv) { time_t systime = time(nullptr); restRun->SetStartTimeStamp((Double_t)systime); - cout << "Events : " << N_events << endl; - if (N_events > 0) // batch mode + cout << "Number of events : " << nEvents << endl; + if (nEvents > 0) // batch mode { G4String command = "/tracking/verbose 0"; UI->ApplyCommand(command); @@ -273,7 +277,7 @@ int main(int argc, char** argv) { UI->ApplyCommand(command); char tmp[256]; - sprintf(tmp, "/run/beamOn %d", N_events); + sprintf(tmp, "/run/beamOn %d", nEvents); command = tmp; UI->ApplyCommand(command); @@ -343,7 +347,7 @@ int main(int argc, char** argv) { } } - else if (N_events == 0) // define visualization and UI terminal for interactive mode + else if (nEvents == 0) // define visualization and UI terminal for interactive mode { cout << "Entering vis mode.." << endl; #ifdef G4UI_USE @@ -357,11 +361,11 @@ int main(int argc, char** argv) { #endif } - else // N_events == -1 + else // nEvents == -1 { - cout << "++++++++++ ERRORRRR +++++++++" << endl; - cout << "++++++++++ ERRORRRR +++++++++" << endl; - cout << "++++++++++ ERRORRRR +++++++++" << endl; + cout << "++++++++++ ERROR +++++++++" << endl; + cout << "++++++++++ ERROR +++++++++" << endl; + cout << "++++++++++ ERROR +++++++++" << endl; cout << "The number of events to be simulated was not recognized properly!" << endl; cout << "Make sure you did not forget the number of events entry in TRestGeant4Metadata." << endl; cout << endl; @@ -369,10 +373,10 @@ int main(int argc, char** argv) { cout << endl; cout << "It should be something like : " << endl; cout << endl; - cout << " " << endl; - cout << "++++++++++ ERRORRRR +++++++++" << endl; - cout << "++++++++++ ERRORRRR +++++++++" << endl; - cout << "++++++++++ ERRORRRR +++++++++" << endl; + cout << " " << endl; + cout << "++++++++++ ERROR +++++++++" << endl; + cout << "++++++++++ ERROR +++++++++" << endl; + cout << "++++++++++ ERROR +++++++++" << endl; cout << endl; } restRun->GetOutputFile()->cd(); @@ -410,7 +414,8 @@ int main(int argc, char** argv) { // writing the geometry object freopen("/dev/null", "w", stdout); freopen("/dev/null", "w", stderr); - REST_Display_CompatibilityMode = true; + + REST_Display_CompatibilityMode = true; // We wait the father process ends properly sleep(5); diff --git a/src/EventAction.cxx b/src/EventAction.cxx index 508b1f86..a7bda9a3 100644 --- a/src/EventAction.cxx +++ b/src/EventAction.cxx @@ -4,6 +4,7 @@ #include #include +#include #include #include #include @@ -16,30 +17,40 @@ extern TRestGeant4Track* restTrack; using namespace std; -EventAction::EventAction() : G4UserEventAction() { restG4Metadata->isFullChainActivated(); } +EventAction::EventAction() : G4UserEventAction() { + fTimer.Start(); + restG4Metadata->isFullChainActivated(); +} EventAction::~EventAction() {} -void EventAction::BeginOfEventAction(const G4Event* geant4_event) { - G4int event_number = geant4_event->GetEventID(); - - restG4Metadata->GetVerboseLevel(); +void EventAction::BeginOfEventAction(const G4Event* event) { + const auto eventID = event->GetEventID(); if (restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Debug) { - cout << "DEBUG: Start of event ID " << event_number << " (" << event_number + 1 << " of " - << restG4Metadata->GetNumberOfEvents() << "). " << restRun->GetEntries() << " Events stored." - << endl; - } else if ((restG4Metadata->PrintProgress() || restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Info) && - geant4_event->GetEventID() % 10000 == 0) { - cout << "INFO: Start of event ID " << event_number << " (" << event_number + 1 << " of " - << restG4Metadata->GetNumberOfEvents() << "). " << restRun->GetEntries() << " Events stored." - << endl - << endl; + G4cout << "DEBUG: Start of event ID " << eventID << " (" << eventID + 1 << " of " + << G4RunManager::GetRunManager()->GetNumberOfEventsToBeProcessed() << "). " + << restRun->GetEntries() << " Events stored." << endl; + } else { + const int numberOfEventsToBePercent = + G4RunManager::GetRunManager()->GetNumberOfEventsToBeProcessed() / 100; + if ((restG4Metadata->PrintProgress() || + restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Essential) && + // print roughly every 1% of events or whenever 10 seconds without printing have elapsed + ((numberOfEventsToBePercent > 0 && (eventID + 1) % numberOfEventsToBePercent == 0) || + fTimer.RealTime() > 10.0)) { + fTimer.Start(); + G4cout << "ESSENTIAL: Start of event ID " << eventID << " (" << eventID + 1 << " of " + << restG4Metadata->GetNumberOfEvents() << "). " << restRun->GetEntries() + << " Events stored" << endl; + } else { + fTimer.Continue(); + } } restTrack->Initialize(); - restG4Event->SetID(event_number); + restG4Event->SetID(eventID); restG4Event->SetOK(true); time_t system_time = time(nullptr); @@ -47,41 +58,45 @@ void EventAction::BeginOfEventAction(const G4Event* geant4_event) { // Defining if the hits in a given volume will be stored for (int i = 0; i < restG4Metadata->GetNumberOfActiveVolumes(); i++) { - Double_t rndNumber = G4UniformRand(); - - if (restG4Metadata->GetStorageChance(i) >= rndNumber) + if (restG4Metadata->GetStorageChance(i) >= 1.00) { restG4Event->ActivateVolumeForStorage(i); - else - restG4Event->DisableVolumeForStorage(i); + } else { + Double_t randomNumber = G4UniformRand(); + if (restG4Metadata->GetStorageChance(i) >= randomNumber) { + restG4Event->ActivateVolumeForStorage(i); + } else { + restG4Event->DisableVolumeForStorage(i); + } + } } } -void EventAction::EndOfEventAction(const G4Event* geant4_event) { - G4int event_number = geant4_event->GetEventID(); +void EventAction::EndOfEventAction(const G4Event* event) { + G4int eventID = event->GetEventID(); if (restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Extreme) { restG4Event->PrintEvent(); } - Double_t minimum_energy_stored = restG4Metadata->GetMinimumEnergyStored(); - Double_t maximum_energy_stored = restG4Metadata->GetMaximumEnergyStored(); + Double_t minimumEnergyStored = restG4Metadata->GetMinimumEnergyStored(); + Double_t maximumEnergyStored = restG4Metadata->GetMaximumEnergyStored(); - int NSubs = SetTrackSubEventIDs(); + int nSubs = SetTrackSubEventIDs(); - Bool_t is_sensitive = false; + Bool_t isSensitive = false; - for (int subId = 0; subId < NSubs + 1; subId++) { + for (int subId = 0; subId < nSubs + 1; subId++) { FillSubEvent(subId); Double_t total_deposited_energy = subRestG4Event->GetTotalDepositedEnergy(); Double_t sensitive_volume_deposited_energy = subRestG4Event->GetSensitiveVolumeEnergy(); - if (minimum_energy_stored < 0) minimum_energy_stored = 0; - if (maximum_energy_stored == 0) maximum_energy_stored = total_deposited_energy + 1.; + if (minimumEnergyStored < 0) minimumEnergyStored = 0; + if (maximumEnergyStored == 0) maximumEnergyStored = total_deposited_energy + 1.; - is_sensitive = - (sensitive_volume_deposited_energy > 0 && total_deposited_energy > minimum_energy_stored && - total_deposited_energy < maximum_energy_stored) || + isSensitive = + (sensitive_volume_deposited_energy > 0 && total_deposited_energy > minimumEnergyStored && + total_deposited_energy < maximumEnergyStored) || restG4Metadata->GetSaveAllEvents(); if (restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Info) { @@ -90,21 +105,22 @@ void EventAction::EndOfEventAction(const G4Event* geant4_event) { debug_level = "DEBUG"; } - if (is_sensitive || restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Debug) { - cout << debug_level - << ": Energy deposited in ACTIVE and SENSITIVE volumes: " << total_deposited_energy - << " keV" << endl; - cout << debug_level - << ": Energy deposited in SENSITIVE volume: " << sensitive_volume_deposited_energy - << " keV" << endl; + if (isSensitive || + restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Debug) { + G4cout << debug_level + << ": Energy deposited in ACTIVE and SENSITIVE volumes: " << total_deposited_energy + << " keV" << endl; + G4cout << debug_level + << ": Energy deposited in SENSITIVE volume: " << sensitive_volume_deposited_energy + << " keV" << endl; } } - if (is_sensitive) { + if (isSensitive) { sensitive_volume_hits_count += 1; // call `ReOrderTrackIds` which before was integrated into `FillSubEvent` - // it takes a while to run so we only do it if we are going to save the event + // it takes a while to run, so we only do it if we are going to save the event ReOrderTrackIds(subId); // fill analysis tree @@ -114,37 +130,40 @@ void EventAction::EndOfEventAction(const G4Event* geant4_event) { analysis_tree->Fill(); } else { // analysis tree is not found (nullptr) - if (restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Warning) { - cout << "WARNING: Analysis tree is not found ('nullptr'). Cannot write event info" - << endl; + if (restG4Metadata->GetVerboseLevel() >= + TRestStringOutput::REST_Verbose_Level::REST_Warning) { + G4cout << "WARNING: Analysis tree is not found ('nullptr'). Cannot write event info" + << endl; } } // fill event tree - TTree* event_tree = restRun->GetEventTree(); - if (event_tree) { - event_tree->Fill(); + TTree* eventTree = restRun->GetEventTree(); + if (eventTree) { + eventTree->Fill(); } else { // event tree is not found (nullptr) - if (restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Warning) { - cout << "WARNING: Event tree is not found ('nullptr'). Cannot write event info" << endl; + if (restG4Metadata->GetVerboseLevel() >= + TRestStringOutput::REST_Verbose_Level::REST_Warning) { + G4cout << "WARNING: Event tree is not found ('nullptr'). Cannot write event info" << endl; } } } } if (restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Info) { - string debug_level = "INFO"; + string debugLevel = "INFO"; if (restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Debug) { - debug_level = "DEBUG"; + debugLevel = "DEBUG"; } - if (is_sensitive || restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Debug) { - cout << debug_level - << ": Events depositing energy in sensitive volume: " << sensitive_volume_hits_count << "/" - << event_number + 1 << endl; - cout << debug_level << ": End of event ID " << event_number << " (" << event_number + 1 << " of " - << restG4Metadata->GetNumberOfEvents() << ")" << endl; - cout << endl; + if (isSensitive || + restG4Metadata->GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Debug) { + G4cout << debugLevel + << ": Events depositing energy in sensitive volume: " << sensitive_volume_hits_count << "/" + << eventID + 1 << endl; + G4cout << debugLevel << ": End of event ID " << eventID << " (" << eventID + 1 << " of " + << G4RunManager::GetRunManager()->GetNumberOfEventsToBeProcessed() << ")" << endl; + G4cout << endl; } } } @@ -183,17 +202,16 @@ void EventAction::FillSubEvent(Int_t subId) { } for (int n = 0; n < restG4Event->GetNumberOfTracks(); n++) { - const auto& tck = restG4Event->GetTrack(n); - - if (tck.GetSubEventID() == subId) - if (tck.GetNumberOfHits() > 0 || restG4Metadata->RegisterEmptyTracks()) { - subRestG4Event->AddTrack(tck); + const auto& track = restG4Event->GetTrack(n); + if (track.GetSubEventID() == subId) { + if (track.GetNumberOfHits() > 0 || restG4Metadata->RegisterEmptyTracks()) { + subRestG4Event->AddTrack(track); } + } } if (restG4Metadata->isVolumeStored(restG4Metadata->GetSensitiveVolume())) { Int_t sensVolID = restG4Metadata->GetActiveVolumeID(restG4Metadata->GetSensitiveVolume()); - subRestG4Event->SetSensitiveVolumeEnergy(subRestG4Event->GetEnergyDepositedInVolume(sensVolID)); } } @@ -206,10 +224,12 @@ void EventAction::ReOrderTrackIds(Int_t subId) { if (subId > 0) { for (int n = 0; n < restG4Event->GetNumberOfTracks(); n++) { - const auto& tck = restG4Event->GetTrack(n); - - if (tck.GetSubEventID() == subId - 1) - if (tck.isRadiactiveDecay()) subRestG4Event->SetSubEventTag(tck.GetParticleName()); + const auto& track = restG4Event->GetTrack(n); + if (track.GetSubEventID() == subId - 1) { + if (track.isRadiactiveDecay()) { + subRestG4Event->SetSubEventTag(track.GetParticleName()); + } + } } } @@ -218,31 +238,36 @@ void EventAction::ReOrderTrackIds(Int_t subId) { Int_t nTracks = subRestG4Event->GetNumberOfTracks(); for (int i = 0; i < nTracks; i++) { - const auto& tr = subRestG4Event->GetTrackPointer(i); - tr->SetTrackID(tr->GetTrackID() - lowestID + 1); - tr->SetParentID(tr->GetParentID() - lowestID + 1); - if (tr->GetParentID() < 0) tr->SetParentID(0); + const auto& track = subRestG4Event->GetTrackPointer(i); + track->SetTrackID(track->GetTrackID() - lowestID + 1); + track->SetParentID(track->GetParentID() - lowestID + 1); + if (track->GetParentID() < 0) { + track->SetParentID(0); + } } for (int i = 0; i < nTracks; i++) { - TRestGeant4Track* tr = subRestG4Event->GetTrackPointer(i); - Int_t id = tr->GetTrackID(); + TRestGeant4Track* track = subRestG4Event->GetTrackPointer(i); + Int_t id = track->GetTrackID(); if (id - i != 1) { // Changing track ids - tr->SetTrackID(i + 1); + track->SetTrackID(i + 1); for (int t = i + 1; t < subRestG4Event->GetNumberOfTracks(); t++) { - TRestGeant4Track* tr2 = subRestG4Event->GetTrackPointer(t); - if (tr2->GetTrackID() == i + 1) tr2->SetTrackID(id); + TRestGeant4Track* trackAux = subRestG4Event->GetTrackPointer(t); + if (trackAux->GetTrackID() == i + 1) { + trackAux->SetTrackID(id); + } } // Changing parent ids for (int t = 0; t < subRestG4Event->GetNumberOfTracks(); t++) { - TRestGeant4Track* tr2 = subRestG4Event->GetTrackPointer(t); - if (tr2->GetParentID() == id) - tr2->SetParentID(i + 1); - else if (tr2->GetParentID() == i + 1) - tr2->SetParentID(id); + TRestGeant4Track* trackAux = subRestG4Event->GetTrackPointer(t); + if (trackAux->GetParentID() == id) { + trackAux->SetParentID(i + 1); + } else if (trackAux->GetParentID() == i + 1) { + trackAux->SetParentID(id); + } } } } @@ -260,42 +285,41 @@ int EventAction::SetTrackSubEventIDs() { } // scan backwards to the parent tracks and set the track's sub event id - int max_subid = 0; + int maxSubEventID = 0; for (auto iter = tracks.begin(); iter != tracks.end(); iter++) { TRestGeant4Track* track = iter->second; - if (track->GetParentID() == 0) { track->SetSubEventID(0); } else { - double tadd = 0; - int parentid = track->GetParentID(); - TRestGeant4Track* ptrack = tracks[parentid]; - while (1) { - if (ptrack) { - tadd += ptrack->GetTrackTimeLength(); - if (tadd > timeDelay) { - int subid = ptrack->GetSubEventID() + 1; - track->SetSubEventID(subid); - if (max_subid < subid) { - max_subid = subid; + double trackTimeLengthTotal = 0; + int parentID = track->GetParentID(); + TRestGeant4Track* trackAux = tracks[parentID]; + while (true) { + if (trackAux) { + trackTimeLengthTotal += trackAux->GetTrackTimeLength(); + if (trackTimeLengthTotal > timeDelay) { + int subEventIDAux = trackAux->GetSubEventID() + 1; + track->SetSubEventID(subEventIDAux); + if (maxSubEventID < subEventIDAux) { + maxSubEventID = subEventIDAux; } break; } else { - if (ptrack->GetParentID() == 0) { + if (trackAux->GetParentID() == 0) { track->SetSubEventID(0); break; } else { - ptrack = tracks[ptrack->GetParentID()]; + trackAux = tracks[trackAux->GetParentID()]; continue; } } } else { - cout << "error! parent track is null" << endl; + G4cout << "error! parent track is null" << endl; abort(); } } } } - return max_subid; + return maxSubEventID; } diff --git a/src/RunAction.cxx b/src/RunAction.cxx index 689d6349..4761c933 100644 --- a/src/RunAction.cxx +++ b/src/RunAction.cxx @@ -22,77 +22,17 @@ RunAction::RunAction(PrimaryGeneratorAction* gen) : G4UserRunAction(), fPrimary( RunAction::~RunAction() {} void RunAction::BeginOfRunAction(const G4Run*) { + G4cout << "========================== Begin of Run Action ========================" << endl; + G4cout << G4RunManager::GetRunManager()->GetNumberOfEventsToBeProcessed() << " events to be simulated" + << endl; + G4cout << "=======================================================================" << endl; // inform the runManager to save random number seed G4RunManager::GetRunManager()->SetRandomNumberStore(false); } -void RunAction::EndOfRunAction(const G4Run* run) { - G4int nbEvents = run->GetNumberOfEvent(); - if (nbEvents == 0) { - return; - } - - G4ParticleDefinition* particle = fPrimary->GetParticleGun()->GetParticleDefinition(); - G4String partName = particle->GetParticleName(); - // G4double eprimary = fPrimary->GetParticleGun()->GetParticleEnergy(); - - G4cout << "======================== run summary ======================"; - G4cout << "\n" << nbEvents << " Events simulated, " << restRun->GetEntries() << " Events stored\n"; - G4cout << "==========================================================="; - G4cout << G4endl; - - // restG4Metadata->PrintMetadata(); - - /* - - G4int prec = 4, wid = prec + 2; - G4int dfprec = G4cout.precision(prec); - - //particle count - // - G4cout << " Nb of generated particles: \n" << G4endl; - - map::iterator it; - for (it = fParticleCount.begin(); it != fParticleCount.end(); it++) { - G4String name = it->first; - G4int count = it->second; - G4double eMean = fEmean[name]/count; - G4double eMin = fEmin[name], eMax = fEmax[name]; - - G4cout << " " << setw(13) << name << ": " << setw(7) << count - << " Emean = " << setw(wid) << G4BestUnit(eMean, "Energy") - << "\t( " << G4BestUnit(eMin, "Energy") - << " --> " << G4BestUnit(eMax, "Energy") - << ")" << G4endl; - } - - //energy momentum balance - // - - if (fDecayCount > 0) { - G4double Ebmean = fEkinTot[0]/fDecayCount; - G4double Pbmean = fPbalance[0]/fDecayCount; - - G4cout << "\n Ekin Total (Q): mean = " - << setw(wid) << G4BestUnit(Ebmean, "Energy") - << "\t( " << G4BestUnit(fEkinTot[1], "Energy") - << " --> " << G4BestUnit(fEkinTot[2], "Energy") - << ")" << G4endl; - - G4cout << "\n Momentum balance (excluding gamma desexcitation): mean = " - << setw(wid) << G4BestUnit(Pbmean, "Energy") - << "\t( " << G4BestUnit(fPbalance[1], "Energy") - << " --> " << G4BestUnit(fPbalance[2], "Energy") - << ")" << G4endl; - } - - // remove all contents in fParticleCount - // - fParticleCount.clear(); - fEmean.clear(); fEmin.clear(); fEmax.clear(); - - // restore default precision - // - G4cout.precision(dfprec); - */ +void RunAction::EndOfRunAction(const G4Run*) { + G4cout << "============================= Run Summary =============================" << endl; + G4cout << restRun->GetEntries() << " events stored out of " + << G4RunManager::GetRunManager()->GetNumberOfEventsToBeProcessed() << " simulated events" << endl; + G4cout << "=======================================================================" << endl; }