diff --git a/.containers/Dockerfile.in b/.containers/Dockerfile.in
new file mode 100644
index 0000000000000000000000000000000000000000..06567d32b0985cde7c1d20a331ba75bda4ab6522
--- /dev/null
+++ b/.containers/Dockerfile.in
@@ -0,0 +1,22 @@
+FROM eicweb.phy.anl.gov:4567/containers/image_recipes/root_spack:@CONTAINER_TAG@
+
+LABEL maintainer "Sylvester Joosten <sjoosten@anl.gov>"
+
+RUN cd /tmp \
+ && git clone https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/hallac_evio.git \
+ && mkdir hallac_evio/build && cd hallac_evio/build \
+ && cmake ../. && make -j20 && make install \
+ && cd /tmp && rm -rf hallac_evio \
+ && git clone https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/analyzer.git \
+ && mkdir analyzer/build && cd analyzer/build && git pull && git checkout v1.8.3 \
+ && cmake ../. && make -j20 VERBOSE=1 && make install \
+ && cd /tmp && rm -rf analyzer
+
+## 2 layers so we can isolate hcana from the rest of the stack
+
+RUN cd /tmp \
+ && git clone https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/hcana.git \
+ && mkdir hcana/build && cd hcana/build \
+ && git pull \
+ && cmake ../. && make -j20 VERBOSE=1 && make install \
+ && cd /tmp && rm -rf hcana
diff --git a/.containers/Makefile b/.containers/Makefile
new file mode 100644
index 0000000000000000000000000000000000000000..d5e8233fff63af2932d816d30bf0cf9c43862c80
--- /dev/null
+++ b/.containers/Makefile
@@ -0,0 +1,64 @@
+# import config.
+# You can change the default config with `make cnf="config_special.env" build`
+cnf ?= config.env
+include $(cnf)
+# exports variables in config.env as environment variables
+export $(shell sed 's/=.*//' $(cnf))
+
+SHELL = bash
+
+# grep the version from the mix file
+VERSION=$(shell bash version.sh)
+TAG_VERSION=$(VERSION)
+LONG_TAG=$(VERSION)-$(PUBLISH_TAG)
+
+# help will output the help for each task
+# thanks to https://marmelab.com/blog/2016/02/29/auto-documented-makefile.html
+.PHONY: help
+
+help: ## This help.
+ @awk 'BEGIN {FS = ":.*?## "} /^[a-zA-Z_-]+:.*?## / {printf "\033[36m%-30s\033[0m %s\n", $$1, $$2}' $(MAKEFILE_LIST)
+
+.DEFAULT_GOAL := help
+
+# ==========================================================================
+#
+build: ## build the image
+ docker build -t $(APP_NAME):$(LONG_TAG) .
+
+build-nc: ## Build the container without caching (from scratch)
+ docker build --no-cache -t $(APP_NAME):$(LONG_TAG) .
+
+# ==========================================================================
+#
+login: ## Auto login to AWS-ECR unsing aws-cli
+ @docker login -u ${CI_REGISTRY_USER} -p ${CI_REGISTRY_PASSWORD} ${CI_REGISTRY}
+ echo "Login COMPLETE"
+
+# ==========================================================================
+#
+
+publish: login ## Publish a tagged container to ECR
+ @echo 'publish $(PUBLISH_TAG) to $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME)'
+ docker tag $(APP_NAME):$(LONG_TAG) $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME):$(PUBLISH_TAG)
+ docker push $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME):$(PUBLISH_TAG)
+ docker rmi $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME):$(PUBLISH_TAG)
+
+publish-version: login ## Publish the `{version}` taged container to ECR
+ @echo 'publish $(VERSION) to $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME):$(VERSION)'
+ docker tag $(APP_NAME):$(LONG_TAG) $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME):$(VERSION)
+ docker push $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME):$(VERSION)
+ docker rmi $(REG_NAME)/$(GL_REG_GROUP)/$(APP_NAME):$(VERSION)
+
+# ==========================================================================
+# remove container from registry
+unpublish: login
+ @echo 'removing $(PUBLISH_TAG)'
+ curl --request DELETE --header "PRIVATE-TOKEN: $(CI_JOB_TOKEN)" "$(REG_API_URL)/$(PUBLISH_TAG)"
+
+# ==========================================================================
+# cleanup docker registry on system used by runner
+cleanup:
+ @echo 'removing $(REG_NAME):$(LONG_TAG)'
+ docker rmi $(REG_NAME):$(LONG_TAG)
+
diff --git a/.containers/config.env b/.containers/config.env
new file mode 100644
index 0000000000000000000000000000000000000000..79a52b2395efce1644765abfad03b063e2c4aa74
--- /dev/null
+++ b/.containers/config.env
@@ -0,0 +1,17 @@
+# Port to run the container
+
+APP_NAME = hcana
+REPO_NAME = hcana
+
+REG_HOST ?= eicweb.phy.anl.gov
+REG_NAME ?= eicweb.phy.anl.gov:4567
+REG_PORT ?= 4567
+REG_URL ?= https://$(REG_HOST)
+
+GL_GROUP = jlab/hallc/analyzer_software
+GL_REG_GROUP = jlab/hallc/analyzer_software/hcana
+GL_REG_NAME = $(REG_NAME)
+
+PUBLISH_TAG = $(HCANA_TAG)
+
+REG_API_URL = $(REG_URL)/api/v4/projects/$(CI_PROJECT_ID)/registry/repositories/49/tags
diff --git a/containers/docker/Makefile b/.containers/docker/Makefile
similarity index 100%
rename from containers/docker/Makefile
rename to .containers/docker/Makefile
diff --git a/.containers/hcana.def.in b/.containers/hcana.def.in
new file mode 100644
index 0000000000000000000000000000000000000000..12ab7e6377bcc4c41f9911595fe6bd51436efa0c
--- /dev/null
+++ b/.containers/hcana.def.in
@@ -0,0 +1,14 @@
+Bootstrap: docker
+From: eicweb.phy.anl.gov:4567/jlab/hallc/analyzer_software/hcana:@HCANA_TAG@
+
+%help
+ singularity container for hcana development
+
+%labels
+ Maintainer "Whitney Armstrong, Sylvester Joosten"
+ Version v1-dev
+
+%post -c /bin/bash
+ echo " -------------------------------------------------"
+ echo " ===> Image setup complete"
+ echo " -------------------------------------------------"
diff --git a/containers/singularity/Singularity.broadwell b/.containers/singularity/Singularity.broadwell
similarity index 100%
rename from containers/singularity/Singularity.broadwell
rename to .containers/singularity/Singularity.broadwell
diff --git a/.containers/version.sh b/.containers/version.sh
new file mode 100644
index 0000000000000000000000000000000000000000..93ef8dc29dfe97a476861893b492f43f49272866
--- /dev/null
+++ b/.containers/version.sh
@@ -0,0 +1,3 @@
+#!/bin/bash
+
+cat ../VERSION
diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index 0aaf4ed384f2b72d30361ea0bd8f11481cc309b2..e2978a196d1e894661c130e22d987988ca0ce226 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -1,84 +1,113 @@
-image: eicweb.phy.anl.gov:4567/containers/image_recipes/root_base:latest
+variables:
+ CONTAINER_TAG: "1.4.0"
+
+default:
+ image: eicweb.phy.anl.gov:4567/containers/image_recipes/ubuntu_dind:latest
+ tags:
+ - silicon
stages:
+ - config
- build
- - build_docker
- - build_sing_img
- - data_replays
- - data_tests
+ - push
+ - deploy
+ - cleanup
-hcana_build:
- stage: build
- tags:
- - sodium
- script:
- - ls -lrth /usr/local/bin
- - bash bin/compile
-
-hcana_docker:
- stage: build_docker
- only:
- - tags
- tags:
- - eic0 docker
- script:
- - docker login -u $CI_REGISTRY_USER -p $CI_REGISTRY_PASSWORD $CI_REGISTRY
- - cd containers/docker && make release
+workflow:
+ rules:
+ - if: '$CI_PIPELINE_SOURCE == "merge_request_event"'
+ - if: '$CI_COMMIT_BRANCH == "master"'
+ - if: '$CI_COMMIT_TAG'
-hcana_singularity:
- tags:
- - singularity
- stage: build_sing_img
- when: manual
- dependencies:
- - hcana_docker
+env:
+ stage: config
script:
- - /bin/bash .gitlabci/setup.sh
- - mkdir -p build
- - pwd
- - cp containers/singularity/Singularity Singularity.hcana
- - cp Singularity.hcana build/.
- - /bin/bash .gitlabci/build.sh Singularity.hcana
- - cp Singularity.hcana.simg build/.
+ - export HCANA_TAG="latest"
+ - |
+ if [ "x${CI_PIPELINE_SOURCE}" == "xmerge_request_event" ]; then
+ export HCANA_TAG="testing-mr-${CI_MERGE_REQUEST_IID}"
+ fi
+ - echo "CI HCANA_TAG for this pipeline set to ${HCANA_TAG}"
+ - echo "HCANA_TAG=$HCANA_TAG" >> hcana.env
artifacts:
- paths:
- - build/Singularity.hcana
- - build/Singularity.hcana.simg
+ reports:
+ dotenv: hcana.env
-elastic_replay:
- when: manual
- tags:
- - eic0 docker
- stage: data_replays
- dependencies:
- - hcana_singularity
+docker:
+ stage: build
+ needs:
+ - env
script:
- - bash tests/replay_elastic_data.sh
- artifacts:
- paths:
- - ROOTfiles/*
- - build/Singularity.hcana
- - build/Singularity.hcana.simg
-
+ - echo $HCANA_TAG
+ - ./.gitlabci/configure.sh .containers/Dockerfile.in
+ - cd .containers
+ - make build-nc
-elastic_test1:
- when: manual
- tags:
- - eic0 docker
- stage: data_tests
- dependencies:
- - elastic_replay
+## for now only publish the stable versions to the registry.
+## if we change this down the line for some CI-time testing
+## make sure to also enable the cleanup job
+publish:
+ stage: push
+ rules:
+ - if: '$CI_COMMIT_BRANCH == "master"'
+ when: manual
+ - if: '$CI_COMMIT_TAG'
+ when: manual
+ needs:
+ - env
+ - docker
script:
- - bash tests/elastic_test.sh
+ - cd .containers
+ - |
+ if [ "x$CI_COMMIT_TAG" != "x" ]; then
+ make publish publish-version
+ else
+ make publish
+ fi
+ retry:
+ max: 2
+ when:
+ - runner_system_failure
+ - stuck_or_timeout_failure
-elastic_test2:
- when: manual
- tags:
- - eic0 docker
- stage: data_tests
- dependencies:
- - elastic_replay
+.singularity:
+ stage: deploy
+ needs:
+ - env
+ - publish
script:
- - bash tests/elastic_test2.sh
+ - ./.gitlabci/configure.sh .containers/hcana.def.in
+ - mkdir -p build
+ - mv .containers/hcana.def build
+ - singularity build build/hcana.sif build/hcana.def
+ retry:
+ max: 2
+ when:
+ - runner_system_failure
+ - stuck_or_timeout_failure
+singularity:latest:
+ extends: .singularity
+ rules:
+ - if: '$CI_COMMIT_BRANCH == "master"'
+ when: manual
+ - if: '$CI_COMMIT_TAG'
+ when: manual
+ artifacts:
+ expire_in: 90 days
+ paths:
+ - build/hcana.sif
+ - build/hcana.def
+purge_image:
+ stage: cleanup
+ needs:
+ - env
+ - publish
+ rules:
+ - if: '$CI_PIPELINE_SOURCE == "merge_request_event"'
+ when: never
+ script:
+ - cd .containers
+ - make cleanup || true
+ - make unpublish || true
diff --git a/.gitlabci/build.sh b/.gitlabci/build.sh
old mode 100644
new mode 100755
diff --git a/.gitlabci/configure.sh b/.gitlabci/configure.sh
new file mode 100755
index 0000000000000000000000000000000000000000..c59d8a41a734e19c5ba32cc1217f12e281e5a56a
--- /dev/null
+++ b/.gitlabci/configure.sh
@@ -0,0 +1,26 @@
+#!/bin/bash
+
+## Configure a CI file based on a template file (first and only
+## argument to this script).
+
+## Known variables that will be substituted:
+## @CONTAINER_TAG@ - docker tag for the ROOT container
+## @HCANA_TAG@ - output tag for the hcana version
+## @HCANA_BRANCH@ - hcana git branch for build
+
+TEMPLATE_FILE=$1
+OUTPUT_FILE=${TEMPLATE_FILE%.*}
+HCANA_BRANCH="master"
+
+if [ -n "${CI_MERGE_REQUEST_SOURCE_BRANCH_NAME}" ] ; then
+ HCANA_BRANCH=$CI_MERGE_REQUEST_SOURCE_BRANCH_NAME
+fi
+
+echo "Configuring CI file: ${TEMPLATE_FILE}"
+echo "Output will be written to: ${OUTPUT_FILE}"
+
+sed "s/@CONTAINER_TAG@/$CONTAINER_TAG/g" $TEMPLATE_FILE | \
+ sed "s/@HCANA_TAG@/$HCANA_TAG/g" | \
+ sed "s/@HCANA_BRANCH@/$HCANA_BRANCH/g" > ${OUTPUT_FILE}
+
+echo "Done"
diff --git a/.gitlabci/setup.sh b/.gitlabci/setup.sh
deleted file mode 100644
index 5c9bbe517eb1dc83dc33b142507764e10fd555d1..0000000000000000000000000000000000000000
--- a/.gitlabci/setup.sh
+++ /dev/null
@@ -1,34 +0,0 @@
-#!/bin/bash
-
-apt-get update && apt-get install -y wget git \
- build-essential \
- squashfs-tools \
- libtool \
- autotools-dev \
- libarchive-dev \
- automake \
- autoconf \
- uuid-dev \
- python3 \
- libssl-dev
-
-
-sed -i -e 's/^Defaults\tsecure_path.*$//' /etc/sudoers
-
-# Check Python
-echo "Python Version:"
-which python
-python --version
-pip install sregistry[all]
-sregistry version
-
-echo "sregistry Version:"
-
-# Install Singularity
-
-cd /tmp && \
- git clone -b vault/release-2.5 https://www.github.com/sylabs/singularity.git
- cd singularity && \
- ./autogen.sh && \
- ./configure --prefix=/usr/local && \
- make -j8 && make install
diff --git a/.gitlabci/sregistry-gitlab.png b/.gitlabci/sregistry-gitlab.png
deleted file mode 100644
index a14e917a20fb32011329af6bc6abd04012acd0d6..0000000000000000000000000000000000000000
Binary files a/.gitlabci/sregistry-gitlab.png and /dev/null differ
diff --git a/.gitlabci/sregistry-gitlab.xcf b/.gitlabci/sregistry-gitlab.xcf
deleted file mode 100644
index f2742162bfc602f92d3d9c5cd573f29e334e042d..0000000000000000000000000000000000000000
Binary files a/.gitlabci/sregistry-gitlab.xcf and /dev/null differ
diff --git a/VERSION b/VERSION
new file mode 100644
index 0000000000000000000000000000000000000000..88c5fb891dcf1d1647d2b84bac0630cf9570d213
--- /dev/null
+++ b/VERSION
@@ -0,0 +1 @@
+1.4.0
diff --git a/bin/compile b/bin/compile
deleted file mode 100755
index ea6387a783b1d0248b1b727ff622ab78247c79cd..0000000000000000000000000000000000000000
--- a/bin/compile
+++ /dev/null
@@ -1,35 +0,0 @@
-#!/bin/bash
-#set -o nounset
-set -o errexit
-startdir=$(pwd)
-export PYTHONPATH=/usr/local/lib:$PYTHONPATH
-export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
-export PATH=/usr/local/bin:$PATH
-source /usr/local/bin/thisroot.sh
-cd /tmp
-git clone https://github.com/fmtlib/fmt.git && cd fmt
-git checkout 6.2.1 && mkdir /tmp/build && cd /tmp/build
-cmake -DBUILD_SHARED_LIBS=TRUE ../fmt
-make -j20 install
-cd /tmp && rm -r /tmp/build && rm -r /tmp/fmt
-cd /tmp
-git clone https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/hallac_evio.git
-mkdir hallac_evio/build && cd hallac_evio/build
-cmake ../. && make -j20 && make install
-cd /tmp && rm -rf hallac_evio
-cd /tmp
-export PYTHONPATH=/usr/local/lib:$PYTHONPATH
-export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
-export PATH=/usr/local/bin:$PATH
-git clone --depth 1 https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/analyzer.git
-mkdir analyzer/build && cd analyzer/build
-cmake ../. && make -j20 VERBOSE=1 && make install
-cd /tmp && rm -rf analyzer
-export PYTHONPATH=/usr/local/lib:$PYTHONPATH
-export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
-export PATH=/usr/local/bin:$PATH
-source /usr/local/bin/thisroot.sh
-cd $startdir
-mkdir build && cd build
-cmake ../. && make -j20 VERBOSE=1 && make install
-
diff --git a/containers/docker/Dockerfile b/containers/docker/Dockerfile
deleted file mode 100644
index f3a412a9af41a6cb8daffc94fea5cd0826f923d3..0000000000000000000000000000000000000000
--- a/containers/docker/Dockerfile
+++ /dev/null
@@ -1,53 +0,0 @@
-FROM eicweb.phy.anl.gov:4567/containers/image_recipes/root_base:latest
-
-LABEL maintainer "Whitney Armstrong <warmstrong@anl.gov>"
-#
-
-RUN ls -lrth /usr/local/lib/lib*.so \
- && export PYTHONPATH=/usr/local/lib:$PYTHONPATH \
- && export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH \
- && export PATH=/usr/local/bin:$PATH \
- && source /usr/local/bin/thisroot.sh \
- && cd /tmp \
- && git clone https://github.com/fmtlib/fmt.git && cd fmt \
- && git checkout 6.2.1 && mkdir /tmp/build && cd /tmp/build \
- && cmake -DBUILD_SHARED_LIBS=TRUE ../fmt \
- && make -j20 install \
- && cd /tmp && rm -r /tmp/build && rm -r /tmp/fmt \
- && cd /tmp \
- && git clone https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/hallac_evio.git \
- && mkdir hallac_evio/build && cd hallac_evio/build \
- && cmake ../. && make -j20 && make install \
- && cd /tmp && rm -rf hallac_evio \
- && cd /tmp \
- && export PYTHONPATH=/usr/local/lib:$PYTHONPATH \
- && export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH \
- && export PATH=/usr/local/bin:$PATH \
- && git clone https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/analyzer.git \
- && mkdir analyzer/build && cd analyzer/build && git pull && git checkout v1.8.2 \
- && cmake ../. && make -j20 VERBOSE=1 && make install \
- && cd /tmp && rm -rf analyzer \
- && export PYTHONPATH=/usr/local/lib:$PYTHONPATH \
- && export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH \
- && export PATH=/usr/local/bin:$PATH \
- && source /usr/local/bin/thisroot.sh \
- && git clone https://eicweb.phy.anl.gov/jlab/hallc/analyzer_software/hcana.git \
- && mkdir hcana/build && cd hcana/build \
- && git pull \
- && cmake ../. && make -j20 VERBOSE=1 && make install \
- && cd /tmp && rm -rf hcana
-
-
-#-DCMAKE_CXX_FLAGS=" -march=haswell -O3 -mfma -malign-data=cacheline -finline-functions "
-#&& wget -O- https://root.cern.ch/download/root_v6.14.06.source.tar.gz | tar -zxvf - \
-#&& mv root-6.14.06 root_master \
-#RUN which c++ && ls -lrth /usr/bin/c++ && cd /tmp/builds/root_build && make -j38 VERBOSE=1 && make install \
-# && cd /tmp && rm -rf /tmp/root_master && rm -rf /tmp/builds/root_build
-
-#RUN useradd -ms /bin/bash -d /opt/user user
-#USER user
-#WORKDIR /opt/bubble_user
-
-##CMD ["-c" ]
-#ENTRYPOINT ["/bin/bash"]
-
diff --git a/containers/docker/Dockerfile.broadwell b/containers/docker/Dockerfile.broadwell
deleted file mode 100644
index 9185d4528387f17217baf733244bdc6bf2e03a54..0000000000000000000000000000000000000000
--- a/containers/docker/Dockerfile.broadwell
+++ /dev/null
@@ -1,58 +0,0 @@
-# ROOT base
-#
-# A container for the latest root
-#
-FROM whit/image_recipes/ubuntu_base:latest
-LABEL maintainer "Whitney Armstrong <warmstrong@anl.gov>"
-#
-
-RUN cd /tmp \
- && wget http://bitbucket.org/eigen/eigen/get/3.3.4.tar.bz2 \
- && tar -xvf 3.3.4.tar.bz2 \
- && cd eigen-* \
- && mkdir build && cd build \
- && cmake ../. -DCMAKE_CXX_FLAGS=" -march=haswell -O3 -mfma -malign-data=cacheline -finline-functions " \
- && make -j10 > /tmp/eigen_build.log && make install
-
-RUN cd /tmp \
- && git clone --depth=1 https://gitlab.cern.ch/CLHEP/CLHEP.git \
- && mkdir -p builds/clhep_build \
- && cd builds/clhep_build \
- && cmake /tmp/CLHEP/. -DCMAKE_CXX_FLAGS=" -march=haswell -O3 -mfma -malign-data=cacheline -finline-functions "\
- && make -j38 install > /tmp/clhep_build.log \
- && cd /tmp && rm -rf /tmp/CLHEP && rm -rf /tmp/builds/clhep_build
-
-RUN cd /tmp \
-&& git clone https://github.com/VcDevel/Vc.git \
-&& cd Vc \
-&& git submodule update --init \
-&& mkdir build && cd build \
-&& cmake -DCMAKE_INSTALL_PREFIX=/usr/local -DBUILD_TESTING=OFF -DTARGET_ARCHITECTURE=broadwell ../. \
-&& make -j30 > /tmp/vc_build.log \
-&& make install
-
-# Build root from the repo master
-RUN cd /tmp \
- && pwd \
- && git clone --depth=1 https://github.com/root-project/root.git root_master \
- && cd /tmp && mkdir -p builds/root_build \
- && cd builds/root_build \
- && cmake ../../root_master/. -Droot7:BOOL=ON -Dcxx17:BOOL=ON -Dfortran:BOOL=ON \
- -Dgdml:BOOL=ON -Dmathmore:BOOL=ON -Dminuit2:BOOL=ON -Dbuiltin_vdt:BOOL=ON -Dbuiltin_veccore:BOOL=ON \
- -Dvc:BOOL=ON -Dbuiltin_vecgeom:BOOL=ON -Dunuran:BOOL=ON \
- && cd /tmp/builds/root_build && make -j38 > /tmp/root_build.log && make install \
- && cd /tmp && rm -rf /tmp/root_master && rm -rf /tmp/builds/root_build
-
-#-DCMAKE_CXX_FLAGS=" -march=haswell -O3 -mfma -malign-data=cacheline -finline-functions "
-#&& wget -O- https://root.cern.ch/download/root_v6.14.06.source.tar.gz | tar -zxvf - \
-#&& mv root-6.14.06 root_master \
-#RUN which c++ && ls -lrth /usr/bin/c++ && cd /tmp/builds/root_build && make -j38 VERBOSE=1 && make install \
-# && cd /tmp && rm -rf /tmp/root_master && rm -rf /tmp/builds/root_build
-
-#RUN useradd -ms /bin/bash -d /opt/user user
-#USER user
-#WORKDIR /opt/bubble_user
-
-##CMD ["-c" ]
-#ENTRYPOINT ["/bin/bash"]
-
diff --git a/containers/docker/README.md b/containers/docker/README.md
deleted file mode 100644
index d6337961f1b7287f56f0b92aa3c87994f04abd55..0000000000000000000000000000000000000000
--- a/containers/docker/README.md
+++ /dev/null
@@ -1,18 +0,0 @@
-Hall A/C Software
-=================
-
-
-Images for running containers for all aspects of the SANE experimental
-analysis.
-
-The starting point is a pre-built image for the ROOT libraries. (ubuntu + ROOT)
-
-Main software libraries:
-
- - `evio`: Built from https://github.com/whit2333/hallac_evio
- - `analyzer`: Hall A analyzer (podd)from https://github.com/whit2333/analyzer
- - `hcana`: Hall C analyzer from https://github.com/whit2333/hcana
-
-These are all built using the super build project `cool_halls` (https://github.com/whit2333/cool_halls)
-
-
diff --git a/containers/docker/config.env b/containers/docker/config.env
deleted file mode 100644
index 631951825ad49e88b62bd773e249dd89b9c01ac0..0000000000000000000000000000000000000000
--- a/containers/docker/config.env
+++ /dev/null
@@ -1,23 +0,0 @@
-# Port to run the container
-PORT=4000
-
-REG_TOKEN ?= ${CI_IMAGE_BUILD_PAT}
-REG_USER ?= whit
-REG_NAME ?= eicweb.phy.anl.gov:4567
-REG_HOST ?= eicweb.phy.anl.gov:4567
-
-# name of alternate build:
-# Dockerfile.$(ALT_NAME) --> $(APP_NAME)_${ALT_NAME}
-ALT_NAME ?= broadwell
-
-APP_NAME = hcana
-REPO_NAME = hcana
-DH_ORG = hallac
-GL_GROUP = jlab/hallc/analyzer_software
-GL_REG_GROUP = jlab/hallc/analyzer_software
-GL_REG_NAME = hcana
-REPO = hcana
-TAG_VERSION = latest
-
-
-
diff --git a/containers/docker/deploy.env b/containers/docker/deploy.env
deleted file mode 100644
index 909aa72672001a5758a7f8c08619c7aa489defa1..0000000000000000000000000000000000000000
--- a/containers/docker/deploy.env
+++ /dev/null
@@ -1,5 +0,0 @@
-# You have to define the values in {}
-#DOCKER_REPO={account-nr}.dkr.ecr.{region}.amazonaws.com
-## optional aws-cli options
-#AWS_CLI_PROFILE={aws-cli-profile}
-#AWS_CLI_REGION={aws-cli-region}
diff --git a/containers/docker/usage.sh b/containers/docker/usage.sh
deleted file mode 100644
index 0b72b4786a77c190f0247a7f566fc50f93c914d6..0000000000000000000000000000000000000000
--- a/containers/docker/usage.sh
+++ /dev/null
@@ -1,26 +0,0 @@
-# INSTALL
-# - copy the files deploy.env, config.env, version.sh and Makefile to your repo
-# - replace the vars in deploy.env
-# - define the version script
-
-# Build the container
-make build
-
-# Build and publish the container
-make release
-
-# Publish a container to AWS-ECR.
-# This includes the login to the repo
-make publish
-
-# Run the container
-make run
-
-# Build an run the container
-make up
-
-# Stop the running container
-make stop
-
-# Build the container with differnt config and deploy file
-make cnf=another_config.env dpl=another_deploy.env build
\ No newline at end of file
diff --git a/containers/docker/version.sh b/containers/docker/version.sh
deleted file mode 100644
index 504899f94c6af4dcb28d4b7e6017c47614478b23..0000000000000000000000000000000000000000
--- a/containers/docker/version.sh
+++ /dev/null
@@ -1,3 +0,0 @@
-#!/bin/bash
-
-echo "0.1"
diff --git a/containers/singularity/Singularity b/containers/singularity/Singularity
deleted file mode 100644
index e57ff64c3eeaacf29147da31449d3b5ca9ba7ca3..0000000000000000000000000000000000000000
--- a/containers/singularity/Singularity
+++ /dev/null
@@ -1,125 +0,0 @@
-Bootstrap: docker
-From: eicweb.phy.anl.gov:4567/jlab/hallc/analyzer_software/hcana:latest
-
-%help
- Hall A/C container.
- Tools:
- - evio : EVIO DAQ data format
- - analyzer : Hall A analyzer (podd)
- - hcana : Hall C analyzer (hcana)
- - root : root version used for the analyzer
- - rootls, rootbrowse, root_config
-
-%labels
- Maintainer "Whitney Armstrong, Sylvester Joosten"
- Version v1.0
-
-%setup -c /bin/bash
- export SINGULARITY_SHELL=/bin/bash
-
-%environment -c /bin/bash
- export PYTHONPATH=/usr/local/lib:${PYTHONPATH}
- export PATH=/usr/local/bin:${PATH}
- export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
- export ROOT_INCLUDE_PATH=/usr/local/include:/usr/local/include/podd:/usr/local/include/hcana
-
-%post -c /bin/bash
- echo "Hello from post"
- echo "Install additional software here"
- source /usr/local/bin/thisroot.sh
- ## libformat and nlohmann json used heavily in new generation replay scripts
- ## libformat
- #cd /tmp && git clone https://github.com/fmtlib/fmt.git && cd fmt && \
- # git checkout 5.3.0 && mkdir /tmp/build && cd /tmp/build && \
- # cmake -DBUILD_SHARED_LIBS=TRUE ../fmt &&
- # make -j4 install && cd /tmp && rm -r /tmp/build && rm -r /tmp/fmt
- ### json
-
-# =======================
-# global
-# =======================
-
-%runscript
- echo "Launching a shell in the Hall A/C singularity container
- exec bash
-
-
-# =======================
-# root
-# =======================
-%apprun root
- root "$@"
-
-%appenv root
- export PYTHONPATH=/usr/local/lib:${PYTHONPATH}
- export PATH=/usr/local/bin:${PATH}
- export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
- export ROOT_INCLUDE_PATH=/usr/local/include/podd:/usr/local/include/hcana
-
-# =======================
-# analyzer
-# =======================
-%apprun analyzer
- analyzer "$@"
-
-%appenv analyzer
- export PYTHONPATH=/usr/local/lib:${PYTHONPATH}
- export PATH=/usr/local/bin:${PATH}
- export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
- export ROOT_INCLUDE_PATH=/usr/local/include/podd:/usr/local/include/hcana
-
-# =======================
-# hcana
-# =======================
-%apphelp hcana
- Run the Hall-C analyzer with same root-style arguments.
-
-%apprun hcana
- source /usr/local/bin/thisroot.sh
- hcana "$@"
-
-%appenv hcana
- export DB_DIR=DBASE
- export PYTHONPATH=/usr/local/lib:${PYTHONPATH}
- export PATH=/usr/local/bin:${PATH}
- export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
- export ROOTSYS=/usr/local
- export ROOT_INCLUDE_PATH=/usr/local/include
- export ROOT_INCLUDE_PATH=/usr/local/include:/usr/local/include/podd:/usr/local/include/hcana
-
-# =======================
-# root-config
-# =======================
-%apprun root_config
- root-config "$@"
-
-%appenv root_config
- export PYTHONPATH=/usr/local/lib:${PYTHONPATH}
- export PATH=/usr/local/bin:${PATH}
- export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
- export ROOT_INCLUDE_PATH=/usr/local/include/podd:/usr/local/include/hcana
-
-
-# =======================
-# rootbrowse
-# =======================
-%apprun rootbrowse
- rootbrowse "$@"
-
-%appenv rootbrowse
- export PYTHONPATH=/usr/local/lib:${PYTHONPATH}
- export PATH=/usr/local/bin:${PATH}
- export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
- export ROOT_INCLUDE_PATH=/usr/local/include/podd:/usr/local/include/hcana
-
-# =======================
-# rootls
-# =======================
-%apprun rootls
- rootls "$@"
-
-%appenv rootls
- export PYTHONPATH=/usr/local/lib:${PYTHONPATH}
- export PATH=/usr/local/bin:${PATH}
- export LD_LIBRARY_PATH=/usr/local/lib:$LD_LIBRARY_PATH
- export ROOT_INCLUDE_PATH=/usr/local/include/podd:/usr/local/include/hcana
diff --git a/podd b/podd
index ecff1455f2868eee65295d4aa6e772eefb2dd941..355b3ce629f8f47f107e2a0fcfa7ee55b3463b6c 160000
--- a/podd
+++ b/podd
@@ -1 +1 @@
-Subproject commit ecff1455f2868eee65295d4aa6e772eefb2dd941
+Subproject commit 355b3ce629f8f47f107e2a0fcfa7ee55b3463b6c
diff --git a/src/THcHelicityScaler.cxx b/src/THcHelicityScaler.cxx
index 13423a6cc88269cffa2c8c9ec9220d9939e29d19..8b6ff54ec2003422908f742b0d20fcdf9c35c47a 100644
--- a/src/THcHelicityScaler.cxx
+++ b/src/THcHelicityScaler.cxx
@@ -5,29 +5,36 @@
~~~
~~~
- THcHelcityScaler *hhelscaler = new THcHelicityScaler("H","HC helicity scalers");
- // hscaler->SetDebugFile("HHelScaler.txt");
- hhelscaler->SetROC(8); // 5 for HMS defaults to 8 for SHMS
- hhelscaler->SetBankID(0x9801); // Will default to this
- gHaEvtHandlers->Add (hhelscaler);
+ THcHelcityScaler *phelscaler = new THcHelicityScaler("P","HC helicity scalers");
+ phelscaler->SetDebugFile("PHelScaler.txt");
+ phelscaler->SetROC(8); // 5 for HMS defaults to 8 for SHMS
+ phelscaler->SetBankID(9801); // Will default to this
+ gHaEvtHandlers->Add (phelscaler);
~~~
-\author
+\authors: S. A. Wood (saw@jlab.org)
+ C. Yero (cyero@jlab.org)
+
*/
+
+#include <bitset>
#include <cstdio>
#include <cstdlib>
#include <cstring>
-#include <iomanip>
#include <iostream>
#include <iterator>
+#include <map>
#include <sstream>
-#include <unordered_map>
#include "TMath.h"
#include "TNamed.h"
#include "TROOT.h"
#include "TString.h"
+#include "fmt/core.h"
+#include "nlohmann/json.hpp"
//#include "THaEvtTypeHandler.h"
+#include "GenScaler.h"
+#include "Scaler3801.h"
#include "THaCodaData.h"
#include "THaEvData.h"
#include "THaGlobals.h"
@@ -36,13 +43,9 @@
#include "THcHelicityScaler.h"
#include "THcParmList.h"
-
#include "Helper.h"
#include "THaVarList.h"
#include "VarDef.h"
-#include "nlohmann/json.hpp"
-
-using nlohmann::json;
using namespace std;
using namespace Decoder;
@@ -56,16 +59,18 @@ static const UInt_t ICUT = 6;
static const UInt_t MAXCHAN = 32;
static const UInt_t defaultDT = 4;
-// Compute Charge Asymmetries
-static std::unordered_map<std::string, Int_t> bcmindex = {
- {"BCM1", 0}, {"BCM2", 2}, {"Unser", 6}, {"BCM4A", 10}, {"BCM4B", 4}, {"BCM4C", 12}, {"1MHz", 8},
-};
-static const Int_t clockindex = 8;
-static const Double_t clockfreq = 1000000.0;
-
THcHelicityScaler::THcHelicityScaler(const char* name, const char* description)
: hcana::ConfigLogging<THaEvtTypeHandler>(name, description), fBankID(9801),
- fUseFirstEvent(kTRUE), fDelayedType(-1) {
+ fUseFirstEvent(kTRUE), fDelayedType(-1), fBCM_Gain(0), fBCM_Offset(0), fBCM_SatOffset(0),
+ fBCM_SatQuadratic(0), fBCM_delta_charge(0), evcount(0), evcountR(0.0), ifound(0),
+ fNormIdx(-1), fNormSlot(-1), dvars(0), dvarsFirst(0), fScalerTree(0), fOnlyBanks(kFALSE),
+ fClockChan(-1), fLastClock(0) {
+
+ fRocSet.clear();
+ fModuleSet.clear();
+
+ //---------------------------------------------------------------------------
+
fROC = -1;
fNScalerChannels = 32;
@@ -75,136 +80,229 @@ THcHelicityScaler::THcHelicityScaler(const char* name, const char* description)
AddEvtType(5);
AddEvtType(6);
AddEvtType(7);
+ AddEvtType(129);
SetDelayedType(129);
}
THcHelicityScaler::~THcHelicityScaler() {
- for (vector<UInt_t*>::iterator it = fDelayedEvents.begin(); it != fDelayedEvents.end(); ++it)
+
+ if (!TROOT::Initialized()) {
+ delete fScalerTree;
+ }
+ Podd::DeleteContainer(scalers);
+ Podd::DeleteContainer(scalerloc);
+
+ for (auto it = fDelayedEvents.begin(); it != fDelayedEvents.end(); ++it)
delete[] * it;
fDelayedEvents.clear();
}
Int_t THcHelicityScaler::End(THaRunBase* runbase) {
+
// Process any delayed events in order received
- const static char* const here = "THcHelicityScaler::End";
- _param_logger->info("{} Analyzing {} delayed helicity scaler events", here,
- fDelayedEvents.size());
- for (auto it = fDelayedEvents.begin(); it != fDelayedEvents.end(); ++it) {
+ for (std::vector<UInt_t*>::iterator it = fDelayedEvents.begin(); it != fDelayedEvents.end();
+ ++it) {
UInt_t* rdata = *it;
AnalyzeBuffer(rdata);
}
+ evNumberR = -1;
- for (auto it = fDelayedEvents.begin(); it != fDelayedEvents.end(); ++it) {
+ for (vector<UInt_t*>::iterator it = fDelayedEvents.begin(); it != fDelayedEvents.end(); ++it)
delete[] * it;
- }
fDelayedEvents.clear();
- // cout << " -- Helicity Scalers -- " << endl;
+ // Write the helicity variables to scaler tree
+ if (fScalerTree) {
+ fScalerTree->Write();
+ }
+
+ // Compute Scaler Asymmetries
+ /*
+ for(Int_t i=0;i<fNScalerChannels;i++) {
+ if(fScalerSums[i]>0.5) {
+ fScalAsymmetry[i] = (fHScalers[0][i]-fHScalers[1][i])/fScalerSums[i];
+ fScalAsymmetryError[i] = 2*TMath::Sqrt(fHScalers[0][i]*fHScalers[1][i]
+ *fScalerSums[i])
+ /(fScalerSums[i]*fScalerSums[i]);
+ } else {
+ fScalAsymmetry[i] = 0.0;
+ fScalAsymmetryError[i] = 0.0;
+ }
+
+ }
+ */
+
+ //------Compute Time Asymmetries-------
+
+ // C.Y. 12/15/2020 Time Asymmetry / Error Calculations (with scaler_current cut)
+
+ if (fQuartetCount <= 1) {
+ fTimeAsymmetry = -1000.;
+ fTimeAsymmetryError = 0.0;
+ } else {
+ fTimeAsymmetry = fTimeAsymSum / fQuartetCount; // normalize asymmetry to total number of
+ // quartets
+ if (fTimeAsymSum2 >= fQuartetCount * TMath::Power(fTimeAsymmetry, 2)) {
+ fTimeAsymmetryError =
+ TMath::Sqrt((fTimeAsymSum2 - fQuartetCount * TMath::Power(fTimeAsymmetry, 2)) /
+ (fQuartetCount * (fQuartetCount - 1)));
+ } else {
+ fTimeAsymmetryError = 0.0;
+ }
+ }
+
+ //------Compute Scaler Asymmetries-------
+
+ // C.Y. 12/15/2020 Scaler Asymmetry / Error Calculations (with scaler_current cut)
+
for (Int_t i = 0; i < fNScalerChannels; i++) {
- if (fScalerSums[i] > 0.5) {
- fAsymmetry[i] = (fHScalers[0][i] - fHScalers[1][i]) / fScalerSums[i];
- fAsymmetryError[i] = 2 * TMath::Sqrt(fHScalers[0][i] * fHScalers[1][i] * fScalerSums[i]) /
- (fScalerSums[i] * fScalerSums[i]);
+
+ if (fQuartetCount <= 1) {
+ fScalAsymmetry[i] = -1000.;
+ fScalAsymmetryError[i] = 0.0;
} else {
- fAsymmetry[i] = 0.0;
- fAsymmetryError[i] = 0.0;
+ fScalAsymmetry[i] =
+ fScalAsymSum[i] / fQuartetCount; // normalize asymmetry to total number of quartets
+ if (fScalAsymSum2[i] >= fQuartetCount * TMath::Power(fScalAsymmetry[i], 2)) {
+ fScalAsymmetryError[i] =
+ TMath::Sqrt((fScalAsymSum2[i] - fQuartetCount * TMath::Power(fScalAsymmetry[i], 2)) /
+ (fQuartetCount * (fQuartetCount - 1)));
+ } else {
+ fScalAsymmetryError[i] = 0.0;
+ }
}
- // printf("%2d %12.0f %12.0f %12.0f %12.8f\n",i,fScalerSums[i],
- // fHScalers[0][i],fHScalers[1][i],
- // fAsymmetry[i]);
}
- // cout << " ---------------------- " << endl;
+
// json dump of helicity charge info
- json j_helicity;
- int run_number = runbase->GetNumber();
+ const int run_number = runbase->GetNumber();
j_helicity["run_number"] = run_number;
- // Compute Charge Asymmetries
+ //------Compute Charge Asymmetries-------
+
+ // C.Y. 12/14/2020 Charge Asymmetry / Error Calculations (with scaler_current cut)
+
+ // Set the helicity scaler module channels for each BCM
std::map<std::string, Int_t> bcmindex;
- bcmindex["BCM1"] = 0;
- bcmindex["BCM2"] = 2;
- // bcmindex["Unser"] = 6;
- bcmindex["BCM4A"] = 10;
- bcmindex["BCM4B"] = 4;
- bcmindex["BCM4C"] = 12;
- // bcmindex["1MHz"] = 8;
- Int_t clockindex = 8;
- Double_t clockfreq = 1000000.0;
-
- Double_t pclock = fHScalers[0][clockindex];
- Double_t mclock = fHScalers[1][clockindex];
- _param_logger->info("{}: -- Beam Charge Asymmetries -- ", here);
+ bcmindex["BCM1_Hel.scal"] = 0;
+ bcmindex["BCM2_Hel.scal"] = 2;
+ bcmindex["Unser_Hel.scal"] = 6;
+ bcmindex["BCM4A_Hel.scal"] = 10;
+ bcmindex["BCM4B_Hel.scal"] = 4;
+ bcmindex["BCM4C_Hel.scal"] = 12;
+ // bcmindex["1MHz_Hel.scal"] = 8;
+
+ // cout << endl << "---------------------- Beam Charge Asymmetries ---------------------- " <<
+ // endl; cout << " BCM Total Charge Beam ON Beam ON Asymmetry" <<
+ // endl; cout << " Name Charge Asymmetry Charge Asymmetry Error" <<
+ // endl;
+
for (Int_t i = 0; i < fNumBCMs; i++) {
- if (bcmindex.find(fBCM_Name[i]) != bcmindex.end()) {
- Int_t index = bcmindex[fBCM_Name[i]];
- Double_t pcounts = fHScalers[0][index];
- Double_t mcounts = fHScalers[1][index];
- // cout << index << " " << fBCM_Name[i] << " " << pcounts << " " << mcounts
- // << " " << fBCM_Gain[i]
- // << " " << fBCM_Offset[i] << endl;
- Double_t pcharge = (pcounts - (pclock / clockfreq) * fBCM_Offset[i]) / fBCM_Gain[i];
- Double_t mcharge = (mcounts - (mclock / clockfreq) * fBCM_Offset[i]) / fBCM_Gain[i];
- fCharge[i] = pcharge + mcharge;
- if (fCharge[i] > 0.0) {
- fChargeAsymmetry[i] = (pcharge - mcharge) / fCharge[i];
+
+ if (fQuartetCount <= 1) {
+ fChargeAsymmetry[i] = -1000.;
+ fChargeAsymmetryError[i] = 0.0;
+ } else {
+ fChargeAsymmetry[i] =
+ fChargeAsymSum[i] / fQuartetCount; // normalize charge asymmetry to total number of
+ // quartets (as the sum is for every quartet)
+
+ if (fChargeAsymSum2[i] >= fQuartetCount * TMath::Power(fChargeAsymmetry[i], 2)) {
+ fChargeAsymmetryError[i] = TMath::Sqrt(
+ (fChargeAsymSum2[i] - fQuartetCount * TMath::Power(fChargeAsymmetry[i], 2)) /
+ (fQuartetCount * (fQuartetCount - 1)));
} else {
- fChargeAsymmetry[i] = 0.0;
+ fChargeAsymmetryError[i] = 0.0;
}
- _param_logger->info("{:10} {:12.2f} uC {:12.2f} ppm", fBCM_Name[i], fCharge[i],
- fChargeAsymmetry[i] * 1e6);
- j_helicity[fBCM_Name[i]] = {{"charge", fCharge[i]},
- {"charge_asymmetry", fChargeAsymmetry[i]}};
}
- }
+ j_helicity[fBCM_Name[i]] = {{"charge", fChargeSum[i]},
+ {"charge_asymmetry", fChargeAsymmetry[i]},
+ {"charge_asymmetry_error", fChargeAsymmetryError}};
- fTime = (pclock + mclock) / clockfreq;
- if (pclock + mclock > 0) {
- fTimeAsymmetry = (pclock - mclock) / (pclock + mclock);
- } else {
- fTimeAsymmetry = 0.0;
+ // printf("%6s %12.2f %12.8f %12.2f %12.8f %12.8f\n",fBCM_Name[i].c_str(),fCharge[i],
+ // fChargeAsymmetry[i],fChargeSum[i],asy,asyerr);
}
- j_helicity["clock"] = {{"time", fTime}, {"time_asymmetry", fTimeAsymmetry}};
- _param_logger->info("{:10} {:12.2f} uC {:12.8f} ppm", "TIME(s)", fTime, fTimeAsymmetry);
+
+ // Compute +/- helicity Times (no BCM cut)
+ Double_t pclock = fHScalers[0][fClockChan];
+ Double_t mclock = fHScalers[1][fClockChan];
+
+ fTimePlus = pclock / fClockFreq;
+ fTimeMinus = mclock / fClockFreq;
+ fTime = (pclock + mclock) / fClockFreq;
+ // if(pclock+mclock>0) {
+ // fTimeAsymmetry = (pclock-mclock)/(pclock+mclock);
+ //} else {
+ // fTimeAsymmetry = 0.0;
+ //}
+ // printf("TIME(s)%12.2f %12.8f %12.2f\n",fTime, fTimeAsymmetry, fTimeSum);
+ j_helicity["clock"] = {{"time", fTime},
+ {"time_asymmetry", fTimeAsymmetry},
+ {"time_asymmetry_error", fTimeAsymmetryError}};
+
+ // Compute Helicity Trigger Asymmetries (no BCM cut)
if (fNTriggersPlus + fNTriggersMinus > 0) {
fTriggerAsymmetry =
((Double_t)(fNTriggersPlus - fNTriggersMinus)) / (fNTriggersPlus + fNTriggersMinus);
} else {
fTriggerAsymmetry = 0.0;
}
- j_helicity["triggers"] = {{"N_triggers", fTime}, {"trigger_asymmetry", fTriggerAsymmetry}};
- _param_logger->info("{}: -- Beam Charge Asymmetries End --", here);
+ j_helicity["triggers"] = {{"N_triggers", fNTriggersPlus + fNTriggersMinus},
+ {"triggers_asymmetry", fTriggerAsymmetry}};
- {
- std::ofstream json_output_file(fmt::format("hel_scalers_{}.json", run_number));
- json_output_file << std::setw(4) << j_helicity << "\n";
- }
+ // write output to
return 0;
}
Int_t THcHelicityScaler::ReadDatabase(const TDatime& date) {
- char prefix[2];
- prefix[0] = 'g';
- prefix[1] = '\0';
- fNumBCMs = 0;
- string bcm_namelist;
- DBRequest list[] = {
- {"NumBCMs", &fNumBCMs, kInt, 0, 1}, {"BCM_Names", &bcm_namelist, kString}, {0}};
+ char prefix[2];
+ prefix[0] = 'g';
+ prefix[1] = '\0';
+ fNumBCMs = 0;
+ DBRequest list[] = {{"NumBCMs", &fNumBCMs, kInt, 0, 1}, {0}};
gHcParms->LoadParmValues((DBRequest*)&list, prefix);
+
if (fNumBCMs > 0) {
fBCM_Gain.resize(fNumBCMs);
fBCM_Offset.resize(fNumBCMs);
+ fBCM_SatOffset.resize(fNumBCMs);
+ fBCM_SatQuadratic.resize(fNumBCMs);
+ fBCM_delta_charge.resize(fNumBCMs);
+
+ string bcm_namelist;
DBRequest list2[] = {{"BCM_Gain", &fBCM_Gain[0], kDouble, (UInt_t)fNumBCMs},
{"BCM_Offset", &fBCM_Offset[0], kDouble, (UInt_t)fNumBCMs},
+ {"BCM_SatQuadratic", &fBCM_SatQuadratic[0], kDouble, (UInt_t)fNumBCMs, 1},
+ {"BCM_SatOffset", &fBCM_SatOffset[0], kDouble, (UInt_t)fNumBCMs, 1},
+ {"BCM_Names", &bcm_namelist, kString},
+ {"BCM_Current_threshold", &fbcm_Current_Threshold, kDouble, 0, 1},
+ {"BCM_Current_threshold_index", &fbcm_Current_Threshold_Index, kInt, 0, 1},
{0}};
+
+ fbcm_Current_Threshold = 0.0;
+ fbcm_Current_Threshold_Index = 0;
+ for (Int_t i = 0; i < fNumBCMs; i++) {
+ fBCM_SatOffset[i] = 0.;
+ fBCM_SatQuadratic[i] = 0.;
+ }
gHcParms->LoadParmValues((DBRequest*)&list2, prefix);
- fBCM_Name = vsplit(bcm_namelist);
+ vector<string> bcm_names = vsplit(bcm_namelist);
+ for (Int_t i = 0; i < fNumBCMs; i++) {
+ fBCM_Name.push_back(bcm_names[i] + "_Hel.scal");
+ fBCM_delta_charge[i] = 0.;
+ }
}
+ fTotalTime = 0.;
+ fPrevTotalTime = 0.;
+ fDeltaTime = -1.;
+
return kOK;
}
+
void THcHelicityScaler::SetDelayedType(int evtype) {
+
/**
* \brief Delay analysis of this event type to end.
*
@@ -218,6 +316,10 @@ void THcHelicityScaler::SetDelayedType(int evtype) {
Int_t THcHelicityScaler::Analyze(THaEvData* evdata) {
+ // C.Y. | THcScalerEvtHandler::Analyze() uses this flag (which is forced to be 1),
+ // but as to why, it is beyond me. For consistency, I have also used it here.
+ Int_t lfirst = 1;
+
if (!IsMyEvent(evdata->GetEvType()))
return -1;
@@ -227,30 +329,84 @@ Int_t THcHelicityScaler::Analyze(THaEvData* evdata) {
EvDump(evdata);
}
+ if (lfirst && !fScalerTree) {
+
+ lfirst = 0;
+
+ // Assign a name to the helicity scaler tree
+ TString sname1 = "TSHel";
+ TString sname2 = sname1 + fName;
+ TString sname3 = fName + " Scaler Data";
+
+ if (fDebugFile) {
+ *fDebugFile << "\nAnalyze 1st time for fName = " << fName << endl;
+ *fDebugFile << sname2 << " " << sname3 << endl;
+ }
+
+ // Create Scaler Tree
+ fScalerTree = new TTree(sname2.Data(), sname3.Data());
+ fScalerTree->SetAutoSave(200000000);
+
+ TString name, tinfo;
+
+ name = "evcount";
+ tinfo = name + "/D";
+ fScalerTree->Branch(name.Data(), &evcountR, tinfo.Data(), 4000);
+
+ name = "evNumber";
+ tinfo = name + "/D";
+ fScalerTree->Branch(name.Data(), &evNumberR, tinfo.Data(), 4000);
+
+ // C.Y. 12/02/2020 Added actual helicity to be stored in scaler tree
+ name = "actualHelicity";
+ tinfo = name + "/D";
+ fScalerTree->Branch(name.Data(), &actualHelicityR, tinfo.Data(), 4000);
+
+ // C.Y. 12/09/2020 Added quartetphase to be stored in scaler tree
+ name = "quartetPhase";
+ tinfo = name + "/D";
+ fScalerTree->Branch(name.Data(), &quartetphaseR, tinfo.Data(), 4000);
+
+ for (size_t i = 0; i < scalerloc.size(); i++) {
+ name = scalerloc[i]->name;
+ tinfo = name + "/D";
+ fScalerTree->Branch(name.Data(), &dvars[i], tinfo.Data(), 4000);
+ }
+ }
+
UInt_t* rdata = (UInt_t*)evdata->GetRawDataBuffer();
if (evdata->GetEvType() == fDelayedType) { // Save this event for processing later
- Int_t evlen = evdata->GetEvLength();
+ Int_t evlen = evdata->GetEvLength();
+
UInt_t* datacopy = new UInt_t[evlen];
fDelayedEvents.push_back(datacopy);
memcpy(datacopy, rdata, evlen * sizeof(UInt_t));
return 1;
- } else { // A normal event
+ }
+
+ else { // A normal event
if (fDebugFile)
*fDebugFile << "\n\nTHcHelicityScaler :: Debugging event type " << dec << evdata->GetEvType()
<< " event num = " << evdata->GetEvNum() << endl
<< endl;
- evNumber = evdata->GetEvNum();
+ evNumber = evdata->GetEvNum();
+ evNumberR = evNumber;
+
Int_t ret;
if ((ret = AnalyzeBuffer(rdata))) {
- //
+ if (fDebugFile)
+ *fDebugFile << "scaler tree ptr 1 " << fScalerTree << endl;
+ if (fDebugFile)
+ *fDebugFile << "ret = " << ret << endl;
}
+
return ret;
}
}
Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
- fNTrigsInBuf = 0;
- const static char* const here = "THcHelicityScaler::AnalyzeBuffer";
+
+ fNTrigsInBuf = 0;
// Parse the data, load local data arrays.
UInt_t* p = (UInt_t*)rdata;
@@ -260,6 +416,7 @@ Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
Int_t evlen = *p + 1;
Int_t ifound = 0;
+
while (p < plast) {
Int_t banklen = *p;
p++; // point to header
@@ -273,8 +430,8 @@ Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
if (fDebugFile)
*fDebugFile << "ROC: " << roc << " " << evlen << " " << *(p - 1) << hex << " " << *p
<< dec << endl;
- // cout << "ROC: " << roc << " " << evlen << " " << *(p-1) << hex << " " <<
- // *p << dec << endl;
+ // cout << "ROC: " << roc << " " << evlen << " " << *(p-1) << hex << " " << *p << dec <<
+ // endl;
if (roc != fROC) { // Not a ROC with helicity scaler
p += *(p - 1) - 1; // Skip to next ROC
}
@@ -289,9 +446,9 @@ Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
// At any point in the bank where the word is not a matching
// header, we stop.
UInt_t tag = (*p >> 16) & 0xffff; // Bank ID (ROC #)
- // UInt_t num = (*p) & 0xff;
- UInt_t* pnext = p + banklen; // Next bank
- p++; // First data word
+ // UInt_t num = (*p) & 0xff;
+ UInt_t* pnext = p + banklen; // Next bank
+ p++; // First data word
// If the bank is not a helicity scaler bank
// or it is not one of the ROC containing helcity scaler data
// skip to the next bank
@@ -299,7 +456,7 @@ Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
if (tag != fBankID) {
p = pnext; // Fall through to end of the above else if
- // cout << " Skipping to next bank" << endl;
+ // cout << " Skipping to next bank" << endl;
} else {
// This is a helicity scaler bank
@@ -307,15 +464,17 @@ Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
Int_t nevents = (banklen - 2) / fNScalerChannels;
// cout << "# of helicity events in bank:" << " " << nevents << endl;
if (nevents > 100) {
- _param_logger->error("{}: Error! Beam off for too long.", here);
+ cout << "Error! Beam off for too long" << endl;
}
fNTrigsInBuf = 0;
+
// Save helcitiy and quad info for THcHelicity
for (Int_t iev = 0; iev < nevents; iev++) { // find number of helicity events in each bank
Int_t index = fNScalerChannels * iev + 1;
+
+ // C.Y. 11/26/2020 This methods extracts the raw helicity information
AnalyzeHelicityScaler(p + index);
- // cout << "H: " << evNumber << endl;
}
}
}
@@ -325,6 +484,7 @@ Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
if (fDebugFile) {
*fDebugFile << "Scaler Header: " << hex << *p << dec;
}
+
if (nskip == 0) {
if (fDebugFile) {
*fDebugFile << endl;
@@ -351,16 +511,18 @@ Int_t THcHelicityScaler::AnalyzeBuffer(UInt_t* rdata) {
}
Int_t THcHelicityScaler::AnalyzeHelicityScaler(UInt_t* p) {
- const static char* const here = "THcHelicityScaler::AnalyzeHelicityScaler";
- Int_t hbits = (p[0] >> 30) & 0x3; // quartet and helcity bits in scaler word
- Bool_t isquartet = (hbits & 2) != 0;
- Int_t ispos = hbits & 1;
- Int_t actualhelicity = 0;
- fHelicityHistory[fNTrigsInBuf] = hbits;
+ const static char* const here = "THcHelicityScaler::AnalyzeHelicityScaler";
+
+ Int_t hbits = (p[0] >> 30) & 0x3; // quartet and helcity bits in scaler word
+ Bool_t isquartet = (hbits & 2) != 0;
+ Int_t ispos = hbits & 1;
+ Int_t actualhelicity = 0;
+ fHelicityHistory[fNTrigsInBuf] = hbits;
fNTrigsInBuf++;
fNTriggers++;
Int_t quartetphase = (fNTriggers - fFirstCycle) % 4;
+
if (fFirstCycle >= -10) {
if (quartetphase == 0) {
Int_t predicted = RanBit30(fRingSeed_reported);
@@ -387,8 +549,7 @@ Int_t THcHelicityScaler::AnalyzeHelicityScaler(UInt_t* p) {
if (isquartet) { // Helicity and quartet signal for next set of scalers
fFirstCycle = fNTriggers - 3;
quartetphase = (fNTriggers - fFirstCycle) % 4;
- //// Helicity at start of quartet is same as last of quartet, so we can start filling the
- /// seed
+ // Helicity at start of quartet is same as last of quartet, so we can start filling the seed
fRingSeed_reported = ((fRingSeed_reported << 1) | ispos) & 0x3FFFFFFF;
fNBits++;
if (fNBits == 30) {
@@ -412,6 +573,7 @@ Int_t THcHelicityScaler::AnalyzeHelicityScaler(UInt_t* p) {
actualhelicity = -actualhelicity;
}
}
+ quartetphase = (quartetphase + 1) % 4;
#else
actualhelicity = (fRingSeed_actual & 1) ? +1 : -1;
if (quartetphase == 1 || quartetphase == 2) {
@@ -422,21 +584,442 @@ Int_t THcHelicityScaler::AnalyzeHelicityScaler(UInt_t* p) {
fRingSeed_actual = 0;
}
+ // C.Y. 12/09/2020 : Pass quartet phase to scaler tree varable
+ quartetphaseR = quartetphase;
+
+ // C.Y. 11/26/2020 The block of code below is used to extract the helicity information from each
+ // channel of the helicity scaler onto a variable to be stored in the scaler tree. For each
+ // channel, each helicity state (+, -, or undefined) is stored in a single varibale. Each helicity
+ // state is tagged to the corresponding scaler read via 'actualHelicityR' which is stored below..
+
+ // C.Y. Assign actual helicity value to a variable to be written to tree (the value may be (0,
+ // +hel (+1), or -hel(-1)) Where actualhelicity=0 is NOT the MPS. It is zero when the actual
+ // helicity has not been determined.
+ actualHelicityR = actualhelicity;
+
+ // C.Y. 11/26/2020 Loop over all 32 scaler channels for a specific helicity scaler module (SIS
+ // 3801)
+ for (Int_t i = 0; i < fNScalerChannels; i++) {
+
+ // C.Y. 11/26/2020 the count expression below gets the scaler raw helicity information (+, -, or
+ // MPS helicity states) for the ith channel
+ Int_t count = p[i] & 0xFFFFFF; // Bottom 24 bits equivalent of scalers->Decode()
+ fScalerChan[i] =
+ count; // pass the helicity raw information to each helicity scaler channel array element
+ }
+
+ // C.Y. 11/26/2020 The block of code below is used to get a cumulative sum of +/- helicity used
+ // for calculation of the cumulative beam charge asymmetry and other quantities in the ::End()
+ // method
if (actualhelicity != 0) {
+
+ // index=0 (+ helicity), index=1 (- helicity)
Int_t hindex = (actualhelicity > 0) ? 0 : 1;
+
+ // C.Y. 11/24/2020 increment the counter for either '+' or '-' helicity states
(actualhelicity > 0) ? (fNTriggersPlus++) : (fNTriggersMinus++);
+
+ // C.Y. 11/24/2020 Loop over all 32 scaler channels for a specific helicity scaler module (SIS
+ // 3801)
for (Int_t i = 0; i < fNScalerChannels; i++) {
- Int_t count = p[i] & 0xFFFFFF; // Bottom 24 bits
+
+ Int_t count = p[i] & 0xFFFFFF; // Bottom 24 bits equivalent of scalers->Decode()
+
+ // Increment either the '+' (hindex=0) or '-' (hindex=1) helicity counts for each [i] scaler
+ // channel channel of a given module
fHScalers[hindex][i] += count;
fScalerSums[i] += count;
}
}
+
+ // Set the helicity scaler clock to define the time
+ fDeltaTime =
+ fScalerChan[fClockChan] /
+ fClockFreq; // total clock counts / clock_frequency (1MHz) for a specific scaler read interval
+ fTotalTime = fPrevTotalTime + fDeltaTime; // cumulative scaler time
+
+ if (fDeltaTime == 0) {
+ _param_logger->error("{}: ******************* Severe Warning ****************************",
+ here);
+ _param_logger->error("{}: In THcHelicityScaler have found fDeltaTime is zero !!", here);
+ _param_logger->error("{}: ******************* Alert DAQ experts ***************************",
+ here);
+ if (fDebugFile)
+ *fDebugFile << " In THcHelicityScaler have found fDeltaTime is zero !! " << endl;
+ }
+
+ fPrevTotalTime =
+ fTotalTime; // set the current total time to the previous time for the upcoming read
+
+ // C.Y. Nov 27, 2020 : (below) code to write the helicity raw data to a variable
+ // and map the variable to the scaler location
+
+ Double_t scal_current = 0;
+ // Loop over each scaler variable from the map
+ for (size_t i = 0; i < scalerloc.size(); i++) {
+ size_t ivar = scalerloc[i]->ivar;
+ size_t idx = scalerloc[i]->index;
+ size_t ichan = scalerloc[i]->ichan;
+
+ // ANALYZE 1ST SCALER READ SEPARATE (There is no previous read before this one)
+ if (evcount == 0) {
+
+ // Loop over the scaler variable (ivar), helicity slot (idx), and slot channel (ichan) -->
+ // idx=0 (only one helicity module per spectrometer arm)
+ if ((ivar < scalerloc.size()) && (idx < scalers.size()) && (ichan < MAXCHAN)) {
+
+ // If fUseFirstEvent=kTRUE (do not skip 1st read)
+ if (fUseFirstEvent) {
+
+ // Write scaler counts
+ if (scalerloc[ivar]->ikind == ICOUNT) {
+ dvars[ivar] = fScalerChan[ichan];
+ dvarsFirst[ivar] = 0.0;
+ }
+
+ // Write scaler time
+ if (scalerloc[ivar]->ikind == ITIME) {
+ dvars[ivar] = fDeltaTime;
+ dvarsFirst[ivar] = 0.0;
+ }
+
+ // Write scaler rate
+ if (scalerloc[ivar]->ikind == IRATE) {
+ dvars[ivar] = (fScalerChan[ichan]) / fDeltaTime;
+ dvarsFirst[ivar] = 0.0;
+ }
+
+ // Write either scaler current or scaler charge
+ if (scalerloc[ivar]->ikind == ICURRENT || scalerloc[ivar]->ikind == ICHARGE) {
+
+ // Get BCM index
+ Int_t bcm_ind = -1;
+ for (Int_t itemp = 0; itemp < fNumBCMs; itemp++) {
+ size_t match = string(scalerloc[ivar]->name.Data()).find(string(fBCM_Name[itemp]));
+ if (match != string::npos) {
+ bcm_ind = itemp;
+ }
+ }
+
+ // Calculate and write scaler current
+ if (scalerloc[ivar]->ikind == ICURRENT) {
+
+ // set default to zero
+ dvars[ivar] = 0.;
+ // Check bcm index and calculate current in a temporary placeholder, "cur_temp"
+ if (bcm_ind != -1) {
+ Double_t cur_temp =
+ ((fScalerChan[ichan]) / fDeltaTime - fBCM_Offset[bcm_ind]) / fBCM_Gain[bcm_ind];
+ cur_temp = cur_temp + fBCM_SatQuadratic[bcm_ind] *
+ TMath::Max(cur_temp - fBCM_SatOffset[i], 0.0);
+ // Assign the calculated scaler current to dvars
+ dvars[ivar] = cur_temp;
+ }
+ // Check which bcm index to place a bcm current threshold cut later on. Assign the the
+ // beam current read to "scal_current" for later use.
+ if (bcm_ind == fbcm_Current_Threshold_Index) {
+ scal_current = dvars[ivar];
+ }
+ }
+
+ // Calculate andd write scaler charge
+ if (scalerloc[ivar]->ikind == ICHARGE) {
+ // Check bcm index and calculate current in a temporary placeholder, "cur_temp", to
+ // determine the charge later on.
+ if (bcm_ind != -1) {
+ Double_t cur_temp =
+ ((fScalerChan[ichan]) / fDeltaTime - fBCM_Offset[bcm_ind]) / fBCM_Gain[bcm_ind];
+ cur_temp =
+ cur_temp +
+ fBCM_SatQuadratic[bcm_ind] *
+ TMath::Power(TMath::Max(cur_temp - fBCM_SatOffset[bcm_ind], 0.0), 2.0);
+ // Calculate the charge for this scaler read
+ fBCM_delta_charge[bcm_ind] = fDeltaTime * cur_temp;
+ // Assign the charge to dvars
+ dvars[ivar] = fBCM_delta_charge[bcm_ind];
+ }
+ }
+ }
+ }
+
+ // If user has set fUseFirstEvent=kFALSE (then use dvarsFirst to store the information for
+ // the 1st event, and leave dvars at 0.0) (The same calculations as above, but assign
+ // dvars=0.0, so that it does not count when written to the scaler tree)
+ else { // ifnotusefirstevent
+
+ if (scalerloc[ivar]->ikind == ICOUNT) {
+ dvarsFirst[ivar] = fScalerChan[ichan];
+ dvars[ivar] = 0.0;
+ }
+
+ if (scalerloc[ivar]->ikind == ITIME) {
+ dvarsFirst[ivar] = fDeltaTime;
+ dvars[ivar] = 0.0;
+ }
+
+ if (scalerloc[ivar]->ikind == IRATE) {
+ dvarsFirst[ivar] = (fScalerChan[ichan]) / fDeltaTime;
+ dvars[ivar] = 0.0;
+ }
+
+ if (scalerloc[ivar]->ikind == ICURRENT || scalerloc[ivar]->ikind == ICHARGE) {
+ Int_t bcm_ind = -1;
+ for (Int_t itemp = 0; itemp < fNumBCMs; itemp++) {
+ size_t match = string(scalerloc[ivar]->name.Data()).find(string(fBCM_Name[itemp]));
+ if (match != string::npos) {
+ bcm_ind = itemp;
+ }
+ }
+
+ if (scalerloc[ivar]->ikind == ICURRENT) {
+
+ dvarsFirst[ivar] = 0.0;
+ dvars[ivar] = 0.0;
+
+ if (bcm_ind != -1) {
+ Double_t cur_temp =
+ ((fScalerChan[ichan]) / fDeltaTime - fBCM_Offset[bcm_ind]) / fBCM_Gain[bcm_ind];
+ cur_temp =
+ cur_temp + fBCM_SatQuadratic[bcm_ind] *
+ TMath::Power(TMath::Max(cur_temp - fBCM_SatOffset[i], 0.0), 2.0);
+ dvarsFirst[ivar] = cur_temp;
+ }
+
+ if (bcm_ind == fbcm_Current_Threshold_Index)
+ scal_current = dvarsFirst[ivar];
+ }
+
+ if (scalerloc[ivar]->ikind == ICHARGE) {
+
+ if (bcm_ind != -1) {
+ Double_t cur_temp =
+ ((fScalerChan[ichan]) / fDeltaTime - fBCM_Offset[bcm_ind]) / fBCM_Gain[bcm_ind];
+ cur_temp =
+ cur_temp +
+ fBCM_SatQuadratic[bcm_ind] *
+ TMath::Power(TMath::Max(cur_temp - fBCM_SatOffset[bcm_ind], 0.0), 2.0);
+ fBCM_delta_charge[bcm_ind] = fDeltaTime * cur_temp;
+ dvarsFirst[ivar] = fBCM_delta_charge[bcm_ind];
+ }
+ }
+ }
+ }
+ } else {
+ _param_logger->warn("{}: ERROR: incorrect index {} {} {}", here, ivar, idx, ichan);
+ }
+ }
+
+ // Calculate Scaler Quantities for ALL OTHER SCALER READS (OTHER THAN THE FIRST READ)
+ else { // evcount != 0
+
+ if ((ivar < scalerloc.size()) && (idx < scalers.size()) && (ichan < MAXCHAN)) {
+
+ if (scalerloc[ivar]->ikind == ICOUNT) {
+ dvars[ivar] = fScalerChan[ichan];
+ }
+
+ if (scalerloc[ivar]->ikind == ITIME) {
+ dvars[ivar] = fDeltaTime;
+ }
+
+ if (scalerloc[ivar]->ikind == IRATE) {
+ dvars[ivar] = fScalerChan[ichan] / fDeltaTime;
+ }
+
+ if (scalerloc[ivar]->ikind == ICURRENT || scalerloc[ivar]->ikind == ICHARGE) {
+ Int_t bcm_ind = -1;
+ for (Int_t itemp = 0; itemp < fNumBCMs; itemp++) {
+ size_t match = string(scalerloc[ivar]->name.Data()).find(string(fBCM_Name[itemp]));
+ if (match != string::npos) {
+ bcm_ind = itemp;
+ }
+ }
+ if (scalerloc[ivar]->ikind == ICURRENT) {
+ dvars[ivar] = 0.;
+
+ if (bcm_ind != -1) {
+
+ if (fDeltaTime > 0) {
+ Double_t cur_temp =
+ (fScalerChan[ichan] / fDeltaTime - fBCM_Offset[bcm_ind]) / fBCM_Gain[bcm_ind];
+ cur_temp =
+ cur_temp +
+ fBCM_SatQuadratic[bcm_ind] *
+ TMath::Power(TMath::Max(cur_temp - fBCM_SatOffset[bcm_ind], 0.0), 2.0);
+ dvars[ivar] = cur_temp;
+ }
+ }
+ if (bcm_ind == fbcm_Current_Threshold_Index)
+ scal_current = dvars[ivar];
+ }
+
+ if (scalerloc[ivar]->ikind == ICHARGE) {
+
+ if (bcm_ind != -1) {
+ dvars[ivar] = 0.0;
+
+ if (fDeltaTime > 0) {
+ Double_t cur_temp =
+ (fScalerChan[ichan] / fDeltaTime - fBCM_Offset[bcm_ind]) / fBCM_Gain[bcm_ind];
+ cur_temp =
+ cur_temp +
+ fBCM_SatQuadratic[bcm_ind] *
+ TMath::Power(TMath::Max(cur_temp - fBCM_SatOffset[bcm_ind], 0.0), 2.0);
+ fBCM_delta_charge[bcm_ind] = fDeltaTime * cur_temp;
+ }
+ dvars[ivar] = fBCM_delta_charge[bcm_ind];
+ }
+ }
+ }
+
+ } else {
+ _param_logger->warn("{}: ERROR: incorrect index {} {} {}", here, ivar, idx, ichan);
+ }
+ }
+ }
+
+ // Analyze Scaler Reads ONLY FOR SCAL_CURRENTS > BCM_CURRENT THRESHOLD
+ //(These variables have the name structure: varibaleName_Cut)
+
+ for (size_t i = 0; i < scalerloc.size(); i++) {
+ size_t ivar = scalerloc[i]->ivar;
+ size_t ichan = scalerloc[i]->ichan;
+
+ if (scalerloc[ivar]->ikind == ICUT + ICOUNT) {
+ if (scal_current > fbcm_Current_Threshold) {
+ dvars[ivar] = fScalerChan[ichan];
+ }
+ }
+
+ if (scalerloc[ivar]->ikind == ICUT + ICHARGE) {
+ Int_t bcm_ind = -1;
+ for (Int_t itemp = 0; itemp < fNumBCMs; itemp++) {
+ size_t match = string(scalerloc[ivar]->name.Data()).find(string(fBCM_Name[itemp]));
+ if (match != string::npos) {
+ bcm_ind = itemp;
+ }
+ }
+ if (scal_current > fbcm_Current_Threshold && bcm_ind != -1) {
+ dvars[ivar] = fBCM_delta_charge[bcm_ind];
+ }
+ }
+
+ if (scalerloc[ivar]->ikind == ICUT + ITIME) {
+ if (scal_current > fbcm_Current_Threshold) {
+ dvars[ivar] = fDeltaTime;
+ }
+ }
+ }
+
+ //--------------------------------------
+
+ // C.Y. 12/13/2020 : Calculate the asymmetries / errors at the end of each quartet (S.A. Wood
+ // approach)
+
+ if (actualhelicity != 0) {
+
+ // Reset fHaveCycle to kFALSE at the start of each quartet (e.g. - + + -, + - - +
+ if (quartetphase == 0) {
+ fHaveCycle[0] = fHaveCycle[1] = fHaveCycle[2] = fHaveCycle[3] = kFALSE;
+ }
+
+ // Check if BCM scaler current is above set threshold
+ if (scal_current > fbcm_Current_Threshold &&
+ (quartetphase == 0 || fHaveCycle[max(quartetphase - 1, 0)])) {
+ fHaveCycle[quartetphase] = kTRUE;
+
+ // Loop over each BCM to get the charge for each cycle of a quartet
+ for (Int_t i = 0; i < fNumBCMs; i++) {
+ fChargeCycle[quartetphase][i] = fBCM_delta_charge[i];
+ }
+
+ // Loop over each Scaler Channel to the the counts for each cycle of a quartet
+ for (Int_t i = 0; i < fNScalerChannels; i++) {
+ fScalCycle[quartetphase][i] = fScalerChan[i];
+ }
+
+ // Set the time for each cycle of the quartet
+ fTimeCycle[quartetphase] = fDeltaTime;
+ }
+
+ // Compute asymmetries at the end of each quartet
+ if (quartetphase == 3 && fHaveCycle[3]) {
+
+ // Loop over BCMs
+ for (Int_t i = 0; i < fNumBCMs; i++) {
+
+ // compute charge asymmetry for each quartet at the end of said quartet
+ Double_t asy =
+ actualhelicity *
+ (fChargeCycle[0][i] + fChargeCycle[3][i] - fChargeCycle[1][i] - fChargeCycle[2][i]) /
+ (fChargeCycle[0][i] + fChargeCycle[3][i] + fChargeCycle[1][i] + fChargeCycle[2][i]);
+
+ fChargeSum[i] +=
+ fChargeCycle[0][i] + fChargeCycle[1][i] + fChargeCycle[2][i] + fChargeCycle[3][i];
+
+ // keep track of sums for proper error calculation
+ fChargeAsymSum[i] += asy;
+ fChargeAsymSum2[i] += asy * asy;
+ }
+
+ //-------
+
+ // Loop over Scaler Channels
+ for (Int_t i = 0; i < fNScalerChannels; i++) {
+
+ // compute scaler asymmetry for each quartet at the end of said quartet
+ Double_t asy = actualhelicity *
+ (fScalCycle[0][i] + fScalCycle[3][i] - fScalCycle[1][i] - fScalCycle[2][i]) /
+ (fScalCycle[0][i] + fScalCycle[3][i] + fScalCycle[1][i] + fScalCycle[2][i]);
+
+ fScalSum[i] += fScalCycle[0][i] + fScalCycle[1][i] + fScalCycle[2][i] + fScalCycle[3][i];
+
+ // keep track of sums for proper error calculation
+ fScalAsymSum[i] += asy;
+ fScalAsymSum2[i] += asy * asy;
+ }
+
+ //-------
+
+ // compute time asymmetry for each quartet at the end of said quartet
+ Double_t asy = actualhelicity *
+ (fTimeCycle[0] + fTimeCycle[3] - fTimeCycle[1] - fTimeCycle[2]) /
+ (fTimeCycle[0] + fTimeCycle[3] + fTimeCycle[1] + fTimeCycle[2]);
+
+ // keep track of total time
+ fTimeSum += fTimeCycle[0] + fTimeCycle[1] + fTimeCycle[2] + fTimeCycle[3];
+
+ // keep track of sums for proper error calculation
+ fTimeAsymSum += asy;
+ fTimeAsymSum2 += asy * asy;
+
+ //------
+
+ // keep track of the total number of quartets
+ fQuartetCount++;
+ }
+ }
+
+ //--------------------------------------
+
+ // increment scaler reads
+ evcount = evcount + 1;
+ evcountR = evcount;
+
+ // clear Genscaler scalers
+ for (size_t j = 0; j < scalers.size(); j++)
+ scalers[j]->Clear("");
+
+ // C.Y. 12/02/2020 Fill Scaler Tree Here
+ if (fScalerTree) {
+ fScalerTree->Fill();
+ }
+
return (0);
}
//_____________________________________________________________________________
Int_t THcHelicityScaler::RanBit30(Int_t ranseed) {
-
UInt_t bit7 = (ranseed & 0x00000040) != 0;
UInt_t bit28 = (ranseed & 0x08000000) != 0;
UInt_t bit29 = (ranseed & 0x10000000) != 0;
@@ -451,9 +1034,13 @@ Int_t THcHelicityScaler::RanBit30(Int_t ranseed) {
//_____________________________________________________________________________
THaAnalysisObject::EStatus THcHelicityScaler::Init(const TDatime& date) {
+ const static char* const here = "THcHelicityScaler::Init";
ReadDatabase(date);
+ const int LEN = 200;
+ char cbuf[LEN];
- fStatus = kOK;
+ fStatus = kOK;
+ fNormIdx = -1;
for (vector<UInt_t*>::iterator it = fDelayedEvents.begin(); it != fDelayedEvents.end(); ++it)
delete[] * it;
@@ -469,6 +1056,166 @@ THaAnalysisObject::EStatus THcHelicityScaler::Init(const TDatime& date) {
fROC = 8; // Default to SHMS crate
}
+ // C.Y. 11/26/2020 : Read In and Parse the variables in the helicity scaler map file
+ TString dfile;
+ dfile = fName + "scaler.txt";
+
+ TString sname0 = "Scalevt";
+ TString sname;
+ sname = "hel" + fName + sname0; // This should be: 'helPScalevt' or 'helHScalevt'
+
+ // Open helicity scaler .dat map file
+ FILE* fi = OpenFile(sname.Data(), date);
+ if (!fi) {
+ _param_logger->warn("{}: Cannot find db file for {} scaler event handler", here, fName);
+ return kFileError;
+ }
+
+ size_t minus1 = -1;
+ size_t pos1;
+ string scomment = "#";
+ string svariable = "variable";
+ string smap = "map";
+ vector<string> dbline;
+
+ while (fgets(cbuf, LEN, fi) != NULL) {
+ std::string sin(cbuf);
+ std::string sinput(sin.substr(0, sin.find_first_of("#")));
+ if (fDebugFile)
+ *fDebugFile << "string input " << sinput << endl;
+ dbline = vsplit(sinput);
+ if (dbline.size() > 0) {
+ pos1 = FindNoCase(dbline[0], scomment);
+ if (pos1 != minus1)
+ continue;
+ pos1 = FindNoCase(dbline[0], svariable);
+
+ if (pos1 != minus1 && dbline.size() > 4) {
+ string sdesc = "";
+ for (size_t j = 5; j < dbline.size(); j++)
+ sdesc = sdesc + " " + dbline[j];
+ UInt_t islot = atoi(dbline[1].c_str());
+ UInt_t ichan = atoi(dbline[2].c_str());
+ UInt_t ikind = atoi(dbline[3].c_str());
+ if (fDebugFile)
+ *fDebugFile << "add var " << dbline[1] << " desc = " << sdesc << " islot= " << islot
+ << " " << ichan << " " << ikind << endl;
+ TString tsname(dbline[4].c_str());
+ TString tsdesc(sdesc.c_str());
+ AddVars(tsname, tsdesc, islot, ichan, ikind);
+ // add extra scaler which is cut on the current
+ if (ikind == ICOUNT || ikind == ITIME || ikind == ICHARGE) {
+ tsname = tsname + "Cut";
+ AddVars(tsname, tsdesc, islot, ichan, ICUT + ikind);
+ }
+ }
+
+ pos1 = FindNoCase(dbline[0], smap);
+ if (fDebugFile)
+ *fDebugFile << "map ? " << dbline[0] << " " << smap << " " << pos1 << " "
+ << dbline.size() << endl;
+ if (pos1 != minus1 && dbline.size() > 6) {
+ Int_t imodel, icrate, islot, inorm;
+ UInt_t header, mask;
+ char cdum[20];
+ sscanf(sinput.c_str(), "%s %d %d %d %x %x %d \n", cdum, &imodel, &icrate, &islot, &header,
+ &mask, &inorm);
+ if ((fNormSlot >= 0) && (fNormSlot != inorm))
+ _param_logger->warn("{}: contradictory norm slot {} {}", here, fNormSlot, inorm);
+ fNormSlot =
+ inorm; // slot number used for normalization. This variable is not used but is checked.
+ Int_t clkchan = -1;
+ Double_t clkfreq = 1;
+ if (dbline.size() > 8) {
+ clkchan = atoi(dbline[7].c_str());
+ clkfreq = 1.0 * atoi(dbline[8].c_str());
+ fClockChan = clkchan;
+ fClockFreq = clkfreq;
+ }
+ if (fDebugFile) {
+ *fDebugFile << "map line " << dec << imodel << " " << icrate << " " << islot << endl;
+ *fDebugFile << " header 0x" << hex << header << " 0x" << mask << dec << " " << inorm
+ << " " << clkchan << " " << clkfreq << endl;
+ }
+
+ switch (imodel) {
+ case 3801: // Hall C Helicity Scalers (SIS 3801)
+ scalers.push_back(new Scaler3801(icrate, islot));
+ if (!fOnlyBanks)
+ fRocSet.insert(icrate);
+ fModuleSet.insert(imodel);
+ break;
+ }
+
+ if (scalers.size() > 0) {
+ UInt_t idx = scalers.size() - 1;
+ // Headers must be unique over whole event, not
+ // just within a ROC
+ scalers[idx]->SetHeader(header, mask);
+ // The normalization slot has the clock in it, so we automatically recognize it.
+ // fNormIdx is the index in scaler[] and
+ // fNormSlot is the slot#, checked for consistency
+ if (clkchan >= 0) {
+ scalers[idx]->SetClock(defaultDT, clkchan, clkfreq);
+ _param_logger->info("{}: Setting scaler clock ... channel = {} ... freq = {}", here,
+ clkchan, clkfreq);
+ if (fDebugFile)
+ *fDebugFile << "Setting scaler clock ... channel = " << clkchan
+ << " ... freq = " << clkfreq << endl;
+ fNormIdx = idx;
+ if (islot != fNormSlot)
+ _param_logger->warn("{}: contradictory norm slot {}", here, islot);
+ }
+ }
+ }
+ }
+
+ } // end while loop
+
+ // can't compare UInt_t to Int_t (compiler warning), so do this
+ nscalers = 0;
+ for (size_t i = 0; i < scalers.size(); i++)
+ nscalers++;
+ // need to do LoadNormScaler after scalers created and if fNormIdx found
+ if (fDebugFile)
+ *fDebugFile << "fNormIdx = " << fNormIdx << endl;
+ if ((fNormIdx >= 0) && fNormIdx < nscalers) {
+ for (Int_t i = 0; i < nscalers; i++) {
+ if (i == fNormIdx)
+ continue;
+ scalers[i]->LoadNormScaler(scalers[fNormIdx]);
+ }
+ }
+
+ // Called after AddVars() has been called
+ DefVars();
+
+ // Verify that the slots are not defined twice
+ for (UInt_t i1 = 0; i1 < scalers.size() - 1; i1++) {
+ for (UInt_t i2 = i1 + 1; i2 < scalers.size(); i2++) {
+ if (scalers[i1]->GetSlot() == scalers[i2]->GetSlot())
+ _param_logger->warn("{}: same slot defined twice", here);
+ }
+ }
+ // Identify indices of scalers[] vector to variables.
+ for (UInt_t i = 0; i < scalers.size(); i++) {
+ for (UInt_t j = 0; j < scalerloc.size(); j++) {
+ if (scalerloc[j]->islot == static_cast<UInt_t>(scalers[i]->GetSlot()))
+ scalerloc[j]->index = i;
+ }
+ }
+
+ if (fDebugFile)
+ *fDebugFile << "THcHelicityScaler:: Name of scaler bank " << fName << endl;
+ for (size_t i = 0; i < scalers.size(); i++) {
+ if (fDebugFile) {
+ *fDebugFile << "Scaler # " << i << endl;
+ scalers[i]->SetDebugFile(fDebugFile);
+ scalers[i]->DebugPrint(fDebugFile);
+ }
+ }
+
+ //------Initialize Helicity Variables / Arrays-----
fNTriggers = 0;
fNTrigsInBuf = 0;
fFirstCycle = -100;
@@ -476,97 +1223,212 @@ THaAnalysisObject::EStatus THcHelicityScaler::Init(const TDatime& date) {
fRingSeed_actual = 0;
fNBits = 0;
fNTriggersPlus = fNTriggersMinus = 0;
- fHScalers[0].resize(fNScalerChannels);
- fHScalers[1].resize(fNScalerChannels);
- fScalerSums.resize(fNScalerChannels);
- fAsymmetry.resize(fNScalerChannels);
- fAsymmetryError.resize(fNScalerChannels);
+ fHScalers[0].resize(fNScalerChannels); //+ helicity scaler counts
+ fHScalers[1].resize(fNScalerChannels); //- helicity scaler counts
+ fScalerSums.resize(fNScalerChannels); // total helicity scaler counts (hel=0, excluded)
+ fScalerChan.resize(fNScalerChannels); // C.Y. 11/26/2020 : added array to store helicity
+ // information per channel
+
for (Int_t i = 0; i < fNScalerChannels; i++) {
- fHScalers[0][i] = 0.0;
- fHScalers[1][i] = 0.0;
- fScalerSums[0] = 0.0;
- fAsymmetry[0] = 0.0;
- fAsymmetryError[0] = 0.0;
+ fHScalers[0][i] = 0.0;
+ fHScalers[1][i] = 0.0;
+ fScalerSums[i] = 0.0;
+ fScalerChan[i] = 0.0;
+ }
+
+ fTimePlus = fTimeMinus = 0.0;
+ fTriggerAsymmetry = 0.0;
+
+ //------C.Y.: 12/13/2020 Initialize Variables for quartet-by-quartet asymmetries (include BCM
+ // current threshold cut)--------- (and error calculation)
+
+ for (Int_t i = 0; i < 4; i++) {
+
+ fChargeCycle[i].resize(fNumBCMs);
+ fScalCycle[i].resize(fNScalerChannels);
+
+ fHaveCycle[i] = kFALSE;
+
+ for (Int_t j = 0; j < fNumBCMs; j++) {
+ fChargeCycle[i][j] = 0.0;
+ }
+
+ for (Int_t j = 0; j < fNScalerChannels; j++) {
+ fScalCycle[i][j] = 0.0;
+ }
+
+ fTimeCycle[i] = 0.0;
}
- fCharge.resize(fNumBCMs);
+ // Initialize quartet counter
+ fQuartetCount = 0.0;
+
+ // Initialize variables for time asymmetry calculation
+ fTimeSum = 0.0;
+ fTimeAsymmetry = 0.0;
+ fTimeAsymmetryError = 0.0;
+ fTimeAsymSum = 0.0;
+ fTimeAsymSum2 = 0.0;
+
+ // Initialize variables for charge asymmetry calculation
+ fChargeSum.resize(fNumBCMs);
fChargeAsymmetry.resize(fNumBCMs);
+ fChargeAsymmetryError.resize(fNumBCMs);
+ fChargeAsymSum.resize(fNumBCMs);
+ fChargeAsymSum2.resize(fNumBCMs);
- fTime = fTimeAsymmetry = 0;
- fTriggerAsymmetry = 0.0;
+ for (Int_t i = 0; i < fNumBCMs; i++) {
+ fChargeSum[i] = 0.0;
+ fChargeAsymmetry[i] = 0.0;
+ fChargeAsymSum[i] = 0.0;
+ fChargeAsymSum2[i] = 0.0;
+ }
+ // Initialize variables for scaler asymmetry calculation
+ fScalSum.resize(fNScalerChannels);
+ fScalAsymmetry.resize(fNScalerChannels);
+ fScalAsymmetryError.resize(fNScalerChannels);
+ fScalAsymSum.resize(fNScalerChannels);
+ fScalAsymSum2.resize(fNScalerChannels);
+
+ for (Int_t i = 0; i < fNScalerChannels; i++) {
+ fScalSum[i] = 0.0;
+ fScalAsymmetry[i] = 0.0;
+ fScalAsymSum[i] = 0.0;
+ fScalAsymSum2[i] = 0.0;
+ }
+
+ //------------------------------------------------------------------------------------------
+
+ // Call MakeParms() to define variables to be used in report file
MakeParms();
return kOK;
}
void THcHelicityScaler::MakeParms() {
- /**
- Put Various helicity scaler results in gHcParms so they can be included in results.
- */
+
+ // Put Various helicity scaler results in gHcParms so they can be included in results.
+
+ // no bcm cut required
gHcParms->Define(Form("g%s_hscaler_plus[%d]", fName.Data(), fNScalerChannels),
- "Plus Helcity Scalers", fHScalers[0][0]);
+ "Plus Helcity Scalers", *(fHScalers[0]));
+
gHcParms->Define(Form("g%s_hscaler_minus[%d]", fName.Data(), fNScalerChannels),
- "Minus Helcity Scalers", fHScalers[1][0]);
- gHcParms->Define(Form("g%s_hscaler_sum[%d]", fName.Data(), fNScalerChannels),
- "Helcity Scalers Sum", fScalerSums[0]);
- gHcParms->Define(Form("g%s_hscaler_asy[%d]", fName.Data(), fNScalerChannels),
- "Helicity Scaler Asymmetry[%d]", fAsymmetry[0]);
- gHcParms->Define(Form("g%s_hscaler_asyerr[%d]", fName.Data(), fNScalerChannels),
- "Helicity Scaler Asymmetry Error[%d]", fAsymmetryError[0]);
+ "Minus Helcity Scalers", *(fHScalers[1]));
+
+ // gHcParms->Define(Form("g%s_hscaler_sum[%d]",fName.Data(),fNScalerChannels),
+ // "Helcity Scalers Sum",*fScalerSums);
+
gHcParms->Define(Form("g%s_hscaler_triggers", fName.Data()), "Total Helicity Scaler Triggers",
fNTriggers);
+
gHcParms->Define(Form("g%s_hscaler_triggers_plus", fName.Data()),
"Positive Helicity Scaler Triggers", fNTriggersPlus);
+
gHcParms->Define(Form("g%s_hscaler_triggers_minus", fName.Data()),
"Negative Helicity Scaler Triggers", fNTriggersMinus);
+
+ gHcParms->Define(Form("g%s_hscaler_trigger_asy", fName.Data()), "Helicity Trigger Asymmetry",
+ fTriggerAsymmetry);
+
+ // gHcParms->Define(Form("g%s_hscaler_time_plus",fName.Data()),
+ // "Positive Helicity Time",fTimePlus);
+
+ // gHcParms->Define(Form("g%s_hscaler_time_minus",fName.Data()),
+ // "Negative Helicity Time",fTimeMinus);
+
+ // bcm cut
+ gHcParms->Define(Form("g%s_hscaler_sum[%d]", fName.Data(), fNScalerChannels),
+ "Helcity Scalers Sum", *fScalSum);
+
+ gHcParms->Define(Form("g%s_hscaler_asy[%d]", fName.Data(), fNScalerChannels),
+ "Helicity Scaler Asymmetry[%d]", *fScalAsymmetry);
+
+ gHcParms->Define(Form("g%s_hscaler_asyerr[%d]", fName.Data(), fNScalerChannels),
+ "Helicity Scaler Asymmetry Error[%d]", *fScalAsymmetryError);
+
gHcParms->Define(Form("g%s_hscaler_charge[%d]", fName.Data(), fNumBCMs), "Helicity Gated Charge",
- fCharge[0]);
+ *fChargeSum);
+
gHcParms->Define(Form("g%s_hscaler_charge_asy[%d]", fName.Data(), fNumBCMs),
- "Helicity Gated Charge Asymmetry", fChargeAsymmetry[0]);
- gHcParms->Define(Form("g%s_hscaler_time", fName.Data()), "Helicity Gated Time (sec)", fTime);
+ "Helicity Gated Charge Asymmetry", *fChargeAsymmetry);
+
+ gHcParms->Define(Form("g%s_hscaler_charge_asyerr[%d]", fName.Data(), fNumBCMs),
+ "Helicity Gated Charge Asymmetry Error", *fChargeAsymmetryError);
+
+ gHcParms->Define(Form("g%s_hscaler_time", fName.Data()), "Helicity Gated Time (sec)", fTimeSum);
+
gHcParms->Define(Form("g%s_hscaler_time_asy", fName.Data()), "Helicity Gated Time Asymmetry",
fTimeAsymmetry);
- gHcParms->Define(Form("g%s_hscaler_trigger_asy", fName.Data()), "Helicity Trigger Asymmetry",
- fTriggerAsymmetry);
+
+ gHcParms->Define(Form("g%s_hscaler_time_asyerr", fName.Data()),
+ "Helicity Gated Time Asymmetry Error", fTimeAsymmetryError);
}
-Double_t THcHelicityScaler::GetPlusCharge(const std::string& name) {
- size_t i = 0;
- for (; i < fNumBCMs; ++i) {
- if (name == fBCM_Name[i]) {
- break;
- }
- }
- auto it = bcmindex.find(name);
+//--------------------C.Y. Sep 20, 2020 : Added Utility Functions--------------------
+
+void THcHelicityScaler::AddVars(TString name, TString desc, UInt_t islot, UInt_t ichan,
+ UInt_t ikind) {
+ if (fDebugFile)
+ *fDebugFile << "C.Y. | Calling THcHelicityScaler::AddVars() " << endl;
+ // need to add fName here to make it a unique variable. (Left vs Right HRS, for example)
+ TString name1 = fName + name;
+ TString desc1 = fName + desc;
+ // We don't yet know the correspondence between index of scalers[] and slots.
+ // Will put that in later.
+ scalerloc.push_back(new HCScalerLoc(name1, desc1, 0, islot, ichan, ikind, scalerloc.size()));
+}
- if (it != bcmindex.end() && i < fNumBCMs) {
- auto index = it->second;
- Double_t clock = fHScalers[0][clockindex];
- Double_t counts = fHScalers[0][index];
- return (counts - (clock / clockfreq) * fBCM_Offset[i]) / fBCM_Gain[i];
+void THcHelicityScaler::DefVars() {
+ const char* const here = "THcHelicityScaler::DefVars()";
+ if (fDebugFile)
+ *fDebugFile << "C.Y. | Calling THcHelicityScaler::DefVars() " << endl;
+ // called after AddVars has finished being called.
+ Nvars = scalerloc.size();
+ if (Nvars == 0)
+ return;
+ dvars.resize(Nvars); // dvars is a member of this class
+ dvarsFirst.resize(Nvars); // dvarsFirst is a member of this class
+ dvars = {0.};
+ dvarsFirst = {0.};
+ if (gHaVars) {
+ if (fDebugFile)
+ *fDebugFile << "THcScalerEVtHandler:: Have gHaVars " << gHaVars << endl;
+ } else {
+ _param_logger->warn("{}: No gHaVars ?! Well, that's a problem !!", here);
+ return;
+ }
+ if (fDebugFile)
+ *fDebugFile << "THcHelicityScaler:: scalerloc size " << scalerloc.size() << endl;
+ const Int_t* count = 0;
+ for (size_t i = 0; i < scalerloc.size(); i++) {
+ gHaVars->DefineByType(scalerloc[i]->name.Data(), scalerloc[i]->description.Data(), &dvars[i],
+ kDouble, count);
}
- _param_logger->warn("THcHelicityScaler::GetPlusCharge: Cannot find scaler {}, return 0.", name);
- return 0.;
}
-Double_t THcHelicityScaler::GetMinusCharge(const std::string& name) {
- size_t i = 0;
- for (; i < fNumBCMs; ++i) {
- if (name == fBCM_Name[i]) {
- break;
- }
+size_t THcHelicityScaler::FindNoCase(const string& sdata, const string& skey) {
+ // Find iterator of word "sdata" where "skey" starts. Case insensitive.
+ string sdatalc, skeylc;
+ sdatalc = "";
+ skeylc = "";
+ for (string::const_iterator p = sdata.begin(); p != sdata.end(); ++p) {
+ sdatalc += tolower(*p);
}
- auto it = bcmindex.find(name);
-
- if (it != bcmindex.end() && i < fNumBCMs) {
- auto index = it->second;
- Double_t clock = fHScalers[1][clockindex];
- Double_t counts = fHScalers[1][index];
- return (counts - (clock / clockfreq) * fBCM_Offset[i]) / fBCM_Gain[i];
+ for (string::const_iterator p = skey.begin(); p != skey.end(); ++p) {
+ skeylc += tolower(*p);
}
- _param_logger->warn("THcHelicityScaler::GetMinusCharge: Cannot find scaler {}, return 0.", name);
- return 0.;
+ if (sdatalc.find(skeylc, 0) == string::npos)
+ return -1;
+ return sdatalc.find(skeylc, 0);
+};
+
+void WriteJson(const std::string& path) {
+ std::ofstream ofile{fmt::format("{}/hel_scalers_{}.json", path, j_helicity["run_number"]};
+ ofile << std::setw(4) << j_helicity << "\n";
}
+//---------------------------------------------------------------------------------
+
ClassImp(THcHelicityScaler)
diff --git a/src/THcHelicityScaler.h b/src/THcHelicityScaler.h
index 650eb4bb200a9fdb6c8422998264174395043d51..461c2ce5fa2b4f34bd3e37b2d36047c23f3c9128 100644
--- a/src/THcHelicityScaler.h
+++ b/src/THcHelicityScaler.h
@@ -11,6 +11,7 @@
#include "Decoder.h"
#include "THaEvtTypeHandler.h"
#include "THaRunBase.h"
+#include "THcScalerEvtHandler.h"
#include "TString.h"
#include "TTree.h"
#include "hcana/Logger.h"
@@ -21,6 +22,7 @@
#include <vector>
class THcHelicity;
+class HCScalerLoc;
class THcHelicityScaler : public hcana::ConfigLogging<THaEvtTypeHandler> {
@@ -28,32 +30,43 @@ public:
THcHelicityScaler(const char*, const char*);
virtual ~THcHelicityScaler();
- Int_t Analyze(THaEvData* evdata);
- Int_t AnalyzeBuffer(UInt_t* rdata);
- Int_t AnalyzeHelicityScaler(UInt_t* p);
+ Int_t Analyze(THaEvData* evdata);
+ Int_t AnalyzeBuffer(UInt_t* rdata);
+ Int_t AnalyzeHelicityScaler(UInt_t* p);
virtual EStatus Init(const TDatime& run_time);
- virtual Int_t ReadDatabase(const TDatime& date);
- virtual Int_t End(THaRunBase* r = 0);
-
- virtual void SetUseFirstEvent(Bool_t b = kFALSE) { fUseFirstEvent = b; }
- virtual void SetDelayedType(int evtype);
- virtual void SetROC(Int_t roc) { fROC = roc; }
- virtual void SetBankID(Int_t bankid) { fBankID = bankid; }
- virtual void SetNScalerChannels(Int_t n) { fNScalerChannels = n; }
- virtual Int_t GetNevents() { return fNTrigsInBuf; }
- virtual Int_t GetNcycles() { return fNTriggers; }
- virtual Int_t GetEvNum() { return evNumber; }
- virtual Int_t* GetHelicityHistoryP() { return fHelicityHistory; }
- virtual Int_t GetReportedSeed() { return fRingSeed_reported; }
- virtual Int_t GetReportedActual() { return fRingSeed_actual; }
- virtual Bool_t IsSeedGood() { return fNBits >= 30; }
+ virtual Int_t ReadDatabase(const TDatime& date);
+ virtual Int_t End(THaRunBase* r = 0);
+
+ virtual void SetUseFirstEvent(Bool_t b = kFALSE) { fUseFirstEvent = b; }
+ virtual void SetDelayedType(int evtype);
+ virtual void SetROC(Int_t roc) { fROC = roc; }
+ virtual void SetBankID(Int_t bankid) { fBankID = bankid; }
+ virtual void SetNScalerChannels(Int_t n) { fNScalerChannels = n; }
+ virtual Int_t GetNevents() { return fNTrigsInBuf; }
+ virtual Int_t GetNcycles() { return fNTriggers; }
+ virtual Int_t GetEvNum() { return evNumber; }
+ virtual Int_t* GetHelicityHistoryP() { return fHelicityHistory; }
+ virtual Int_t GetReportedSeed() { return fRingSeed_reported; }
+ virtual Int_t GetReportedActual() { return fRingSeed_actual; }
+ virtual Bool_t IsSeedGood() { return fNBits >= 30; }
Double_t GetPlusCharge(const std::string& name);
Double_t GetMinusCharge(const std::string& name);
+ // utility function to write out json helicity file
+ void WriteJson(const std::string& path);
+
private:
+ //------------C.Y. Sep 20, 2020 :: Added Utility Function Prototypes----------------
+ void AddVars(TString name, TString desc, UInt_t iscal, UInt_t ichan, UInt_t ikind);
+ void DefVars();
+ static size_t FindNoCase(const std::string& sdata, const std::string& skey);
+
+ std::vector<Decoder::GenScaler*> scalers;
+ std::vector<HCScalerLoc*> scalerloc;
+
Int_t RanBit30(Int_t ranseed);
- void MakeParms();
+ void MakeParms();
UInt_t fBankID;
// Helicity Scaler variables
@@ -76,9 +89,44 @@ private:
std::vector<Double_t> fScalerSums, fAsymmetry, fAsymmetryError;
std::vector<Double_t> fCharge, fChargeAsymmetry;
Double_t fTime;
- Double_t fTimeAsymmetry;
+ Double_t fTimePlus;
+ Double_t fTimeMinus;
Double_t fTriggerAsymmetry;
+ //---- C.Y.: 12/14/2020 Variables for quartet-by-quartet asymmetry/error calculations ----
+ Bool_t fHaveCycle[4];
+
+ Int_t fQuartetCount; // keep track of number of quartets
+
+ // quartet-by-quartet time asymmetry variables
+ Double_t fTimeCycle[4];
+ Double_t fTimeSum;
+ Double_t fTimeAsymmetry;
+ Double_t fTimeAsymmetryError;
+ Double_t fTimeAsymSum;
+ Double_t fTimeAsymSum2;
+
+ // quartet-by-quartet scaler counts asymmetry variables
+ std::vector<Double_t> fScalCycle[4];
+ std::vector<Double_t> fScalSum; // reminder: need to initialize
+ std::vector<Double_t> fScalAsymmetry;
+ std::vector<Double_t> fScalAsymmetryError;
+ std::vector<Double_t> fScalAsymSum;
+ std::vector<Double_t> fScalAsymSum2;
+
+ // quartet-by-quartet charge asymmetry variables
+ std::vector<Double_t> fChargeCycle[4];
+ std::vector<Double_t> fChargeSum;
+ std::vector<Double_t> fChargeAsymmetry;
+ std::vector<Double_t> fChargeAsymmetryError;
+ std::vector<Double_t> fChargeAsymSum;
+ std::vector<Double_t> fChargeAsymSum2;
+
+ //----------------------
+
+ //----C.Y. Nov 26, 2020----
+ std::vector<Double_t> fScalerChan;
+
std::vector<UInt_t*> fDelayedEvents;
Int_t fROC;
Int_t fNScalerChannels; // Number of scaler channels/event
@@ -86,12 +134,43 @@ private:
Int_t fNumBCMs;
std::vector<Double_t> fBCM_Gain, fBCM_Offset;
std::vector<std::string> fBCM_Name;
+ //---C.Y. Sep 2020 : Added additional BCM-related variables--
+ std::vector<Double_t> fBCM_SatOffset;
+ std::vector<Double_t> fBCM_SatQuadratic;
+ std::vector<Double_t> fBCM_delta_charge;
+ Double_t fTotalTime;
+ Double_t fDeltaTime;
+ Double_t fPrevTotalTime;
+ Double_t fbcm_Current_Threshold;
+ Double_t fClockFreq;
+ Int_t fbcm_Current_Threshold_Index;
+
+ //----C.Y. Sep 20, 2020 : Added additional variables-----
+ // (required by utility functions and scaler tree output)
+ UInt_t evcount;
+ Double_t evcountR;
+ UInt_t evNumber;
+ Double_t evNumberR;
+ Double_t actualHelicityR;
+ Double_t quartetphaseR;
+ Int_t Nvars, ifound, fNormIdx, fNormSlot, nscalers;
+ std::vector<Double_t> dvars;
+ std::vector<Double_t> dvarsFirst;
+ TTree* fScalerTree;
+ Bool_t fOnlySyncEvents;
+ Bool_t fOnlyBanks;
+ Int_t fClockChan;
+ UInt_t fLastClock;
+ std::set<UInt_t> fRocSet;
+ std::set<UInt_t> fModuleSet;
+ //--------------------------------------------------------
+
+ // json file with helicity info
+ json j_helicity;
THcHelicityScaler(const THcHelicityScaler& fh);
THcHelicityScaler& operator=(const THcHelicityScaler& fh);
- UInt_t evNumber;
-
ClassDef(THcHelicityScaler, 0) // Scaler Event handler
};
diff --git a/src/THcParmList.cxx b/src/THcParmList.cxx
index c2a3c88a47b772a290968e65952403f8c35b186f..deef2e9a9c31ef700685f43810d74ceaf78f40e3 100644
--- a/src/THcParmList.cxx
+++ b/src/THcParmList.cxx
@@ -43,50 +43,48 @@ An instance of THaTextvars is created to hold the string parameters.
#include <iomanip>
#include "TBufferJSON.h"
-#include "nlohmann/json.hpp"
#include "TObjArray.h"
#include "TObjString.h"
#include "TSystem.h"
+#include "nlohmann/json.hpp"
-#include "THcParmList.h"
-#include "THaVar.h"
#include "THaFormula.h"
+#include "THaVar.h"
+#include "THcParmList.h"
#include "TMath.h"
/* #incluce <algorithm> include <fstream> include <cstring> */
-#include <iostream>
-#include <fstream>
#include <cassert>
#include <cstdlib>
-#include <stdexcept>
+#include <fstream>
+#include <iostream>
#include <memory>
+#include <stdexcept>
using namespace std;
-Int_t fDebug = 1; // Keep this at one while we're working on the code
+Int_t fDebug = 1; // Keep this at one while we're working on the code
#if __cplusplus < 201103L
-# define SMART_PTR auto_ptr
+#define SMART_PTR auto_ptr
#else
-# define SMART_PTR unique_ptr
+#define SMART_PTR unique_ptr
#endif
ClassImp(THcParmList)
-/// Create empty numerical and string parameter lists
-THcParmList::THcParmList() : hcana::ConfigLogging<THaVarList>()
-{
+ /// Create empty numerical and string parameter lists
+ THcParmList::THcParmList()
+ : hcana::ConfigLogging<THaVarList>() {
TextList = new THaTextvars;
}
-inline static bool IsComment( const string& s, string::size_type pos )
-{
- return ( pos != string::npos && pos < s.length() &&
- (s[pos] == '#' || s[pos] == ';' || s.substr(pos,2) == "//") );
+inline static bool IsComment(const string& s, string::size_type pos) {
+ return (pos != string::npos && pos < s.length() &&
+ (s[pos] == '#' || s[pos] == ';' || s.substr(pos, 2) == "//"));
}
-void THcParmList::Load( const char* fname, Int_t RunNumber )
-{
+void THcParmList::Load(const char* fname, Int_t RunNumber) {
/**
Most parameter files used in the ENGINE should work.
@@ -127,138 +125,136 @@ The ENGINE CTP support parameter "blocks" which were marked with
static const char* const whtspc = " \t";
- ifstream ifiles[100]; // Should use stack instead
+ ifstream ifiles[100]; // Should use stack instead
- Int_t nfiles=0;
+ Int_t nfiles = 0;
ifiles[nfiles].open(fname);
- if(ifiles[nfiles].is_open()) {
- //cout << "Opening parameter file: [" << nfiles << "] " << fname << endl;
+ if (ifiles[nfiles].is_open()) {
+ // cout << "Opening parameter file: [" << nfiles << "] " << fname << endl;
nfiles++;
}
- if(!nfiles) {
- static const char* const here = "THcParmList::LoadFromFile";
- //Error (here, "error opening parameter file %s",fname);
- return; // Need a success argument returned
+ if (!nfiles) {
+ static const char* const here = "THcParmList::LoadFromFile";
+ // Error (here, "error opening parameter file %s",fname);
+ return; // Need a success argument returned
}
string line;
- char varname[100];
- Int_t InRunRange;
- Int_t currentindex = 0;
- Int_t linecount=0; // Count of non comment/blank lines
+ char varname[100];
+ Int_t InRunRange;
+ Int_t currentindex = 0;
+ Int_t linecount = 0; // Count of non comment/blank lines
varname[0] = '\0';
- if(RunNumber > 0) {
- InRunRange = 0; // Wait until run number range matching RunNumber is found
- //cout << "Reading Parameters for run " << RunNumber << endl;
+ if (RunNumber > 0) {
+ InRunRange = 0; // Wait until run number range matching RunNumber is found
+ // cout << "Reading Parameters for run " << RunNumber << endl;
} else {
- InRunRange = 1; // Interpret all lines
+ InRunRange = 1; // Interpret all lines
}
- while(nfiles) {
+ while (nfiles) {
string current_comment("");
// EJB_Note: existing_comment is never used.
// string existing_comment("");
string::size_type start, pos = 0;
- if(!getline(ifiles[nfiles-1],line)) {
- ifiles[nfiles-1].close();
+ if (!getline(ifiles[nfiles - 1], line)) {
+ ifiles[nfiles - 1].close();
nfiles--;
// cout << nfiles << ": " << "Closed" << endl;
continue;
}
// Look for include statement
- if(line.compare(0,strlen(INCLUDESTR),INCLUDESTR)==0) {
- line.erase(0,strlen(INCLUDESTR));
+ if (line.compare(0, strlen(INCLUDESTR), INCLUDESTR) == 0) {
+ line.erase(0, strlen(INCLUDESTR));
pos = line.find_first_not_of(whtspc);
// Strip leading white space
- if(pos != string::npos && pos > 0 && pos < line.length()) {
- line.erase(0,pos);
+ if (pos != string::npos && pos > 0 && pos < line.length()) {
+ line.erase(0, pos);
}
- char quotechar=line[0];
- if(quotechar == '"' || quotechar == '\'') {
- line.erase(0,1);
- line.erase(line.find_first_of(quotechar));
+ char quotechar = line[0];
+ if (quotechar == '"' || quotechar == '\'') {
+ line.erase(0, 1);
+ line.erase(line.find_first_of(quotechar));
} else {
- line.erase(line.find_first_of(whtspc));
+ line.erase(line.find_first_of(whtspc));
}
// cout << line << endl;
ifiles[nfiles].open(line.c_str());
- if(ifiles[nfiles].is_open()) {
+ if (ifiles[nfiles].is_open()) {
_logger->info("Opening parameter file: [{}] {} ", nfiles, line);
- //cout << "Opening parameter file: [" << nfiles << "] " << line << endl;
- nfiles++;
+ // cout << "Opening parameter file: [" << nfiles << "] " << line << endl;
+ nfiles++;
}
continue;
}
// Blank line or comment?
- if( line.empty()
- || (start = line.find_first_not_of( whtspc )) == string::npos
- || IsComment(line, start) )
+ if (line.empty() || (start = line.find_first_not_of(whtspc)) == string::npos ||
+ IsComment(line, start))
continue;
// Get rid of trailing comments and leading and trailing whitespace
// Need to save the comment and put it in the thVar
- while( (pos = line.find_first_of("#;/", pos+1)) != string::npos ) {
- if( IsComment(line, pos) ) {
- current_comment.assign(line,pos+1,line.length());
- line.erase(pos); // Strip off comment
- // Strip leading white space from comment
- //cout << "CommentA: " << current_comment << endl;
- pos = current_comment.find_first_not_of(whtspc);
- if(pos!=string::npos && pos > 0 && pos < current_comment.length()) {
- current_comment.erase(0,pos);
- }
- //cout << "CommentB: " << current_comment << endl;
- break;
+ while ((pos = line.find_first_of("#;/", pos + 1)) != string::npos) {
+ if (IsComment(line, pos)) {
+ current_comment.assign(line, pos + 1, line.length());
+ line.erase(pos); // Strip off comment
+ // Strip leading white space from comment
+ // cout << "CommentA: " << current_comment << endl;
+ pos = current_comment.find_first_not_of(whtspc);
+ if (pos != string::npos && pos > 0 && pos < current_comment.length()) {
+ current_comment.erase(0, pos);
+ }
+ // cout << "CommentB: " << current_comment << endl;
+ break;
}
}
- pos = line.find_last_not_of( whtspc );
- assert( pos != string::npos );
- if( pos != string::npos && ++pos < line.length() )
+ pos = line.find_last_not_of(whtspc);
+ assert(pos != string::npos);
+ if (pos != string::npos && ++pos < line.length())
line.erase(pos);
pos = line.find_first_not_of(whtspc);
// Strip leading white space
- if(pos != string::npos && pos > 0 && pos < line.length()) {
- line.erase(0,pos);
+ if (pos != string::npos && pos > 0 && pos < line.length()) {
+ line.erase(0, pos);
}
// Ignore begin and end statements
- if(line.compare(0,5,"begin")==0 ||
- line.compare(0,3,"end")==0) {
+ if (line.compare(0, 5, "begin") == 0 || line.compare(0, 3, "end") == 0) {
cout << "Skipping: " << line << endl;
continue;
}
// Get rid of all white space not in quotes
// Step through one char at a time
- pos = 0;
- int inquote=0;
- char quotechar=' ';
+ pos = 0;
+ int inquote = 0;
+ char quotechar = ' ';
// cout << "Unstripped line: |" << line << "|" << endl;
- while(pos<line.length()) {
- if(inquote) {
- if(line[pos++] == quotechar) { // Possibly end of quoted string
- if(line[pos] == quotechar) { // Protected quote
- pos++; // Skip the protected quote
- } else { // End of quoted string
- inquote = 0;
- quotechar = ' ';
- pos++;
- }
- }
+ while (pos < line.length()) {
+ if (inquote) {
+ if (line[pos++] == quotechar) { // Possibly end of quoted string
+ if (line[pos] == quotechar) { // Protected quote
+ pos++; // Skip the protected quote
+ } else { // End of quoted string
+ inquote = 0;
+ quotechar = ' ';
+ pos++;
+ }
+ }
} else {
- if(line[pos] == ' ' || line[pos] == '\t') {
- line.erase(pos,1);
- } else if(line[pos] == '"' || line[pos] == '\'') {
- quotechar = line[pos++];\
- inquote = 1;
- } else {
- pos++;
- }
+ if (line[pos] == ' ' || line[pos] == '\t') {
+ line.erase(pos, 1);
+ } else if (line[pos] == '"' || line[pos] == '\'') {
+ quotechar = line[pos++];
+ inquote = 1;
+ } else {
+ pos++;
+ }
}
}
// cout << "Stripped line: |" << line << "|" << endl;
@@ -269,96 +265,97 @@ The ENGINE CTP support parameter "blocks" which were marked with
linecount++;
// If RunNumber>0 and first line we encounter is not a run range, need to
// print an error
- if(RunNumber>0 && nfiles==1) {
- if(line.find_first_not_of("0123456789-,")==string::npos) { // Interpret as runnum range
- // Interpret line as a list of comma separated run numbers or ranges
- TString runnums(line.c_str());
- SMART_PTR<TObjArray> runnumarr( runnums.Tokenize(",") );
- Int_t nranges=runnumarr->GetLast()+1;
-
- InRunRange = 0;
- Int_t ind;
- for(Int_t i=0;i<nranges;i++) {
- TString runstr = ((TObjString *)runnumarr->At(i))->GetString();
- if(runstr.IsDec()) { // A single run number
- if(RunNumber == runstr.Atoi()) {
- InRunRange = 1;
- break;
- }
- } else if ((ind=runstr.First('-'))>=0) { // A run range
- TString start=runstr(0,ind);
- TString end=runstr(ind+1,runstr.Length());
- if(start.IsDec() && end.IsDec()) {
- if((RunNumber >= start.Atoi()) && (RunNumber <= end.Atoi())) {
- InRunRange = 1;
- break;
- }
- }
- }
- }
- continue; // Skip to next line
+ if (RunNumber > 0 && nfiles == 1) {
+ if (line.find_first_not_of("0123456789-,") == string::npos) { // Interpret as runnum range
+ // Interpret line as a list of comma separated run numbers or ranges
+ TString runnums(line.c_str());
+ SMART_PTR<TObjArray> runnumarr(runnums.Tokenize(","));
+ Int_t nranges = runnumarr->GetLast() + 1;
+
+ InRunRange = 0;
+ Int_t ind;
+ for (Int_t i = 0; i < nranges; i++) {
+ TString runstr = ((TObjString*)runnumarr->At(i))->GetString();
+ if (runstr.IsDec()) { // A single run number
+ if (RunNumber == runstr.Atoi()) {
+ InRunRange = 1;
+ break;
+ }
+ } else if ((ind = runstr.First('-')) >= 0) { // A run range
+ TString start = runstr(0, ind);
+ TString end = runstr(ind + 1, runstr.Length());
+ if (start.IsDec() && end.IsDec()) {
+ if ((RunNumber >= start.Atoi()) && (RunNumber <= end.Atoi())) {
+ InRunRange = 1;
+ break;
+ }
+ }
+ }
+ }
+ continue; // Skip to next line
} else {
- if(linecount==1) {
- _logger->warn("THcParmList::Load in database mode but first line is not\n"
- " a run number or run number range. Parameter definitions\n"
- " will be ignored until a run number or range is specified.\n");
- }
+ if (linecount == 1) {
+ _logger->warn("THcParmList::Load in database mode but first line is not\n"
+ " a run number or run number range. Parameter definitions\n"
+ " will be ignored until a run number or range is specified.\n");
+ }
}
}
- if(!InRunRange) continue;
+ if (!InRunRange)
+ continue;
// Interpret left of = as var name
- Int_t valuestartpos=0; // Stays zero if no = found
- Int_t ttype = 0; // Are any of the values floating point?
- if((pos=line.find_first_of("="))!=string::npos) {
- strcpy(varname, (line.substr(0,pos)).c_str());
- valuestartpos = pos+1;
- currentindex = 0;
+ Int_t valuestartpos = 0; // Stays zero if no = found
+ Int_t ttype = 0; // Are any of the values floating point?
+ if ((pos = line.find_first_of("=")) != string::npos) {
+ strcpy(varname, (line.substr(0, pos)).c_str());
+ valuestartpos = pos + 1;
+ currentindex = 0;
}
// If first char after = is a quote, then this is a string assignment
- if(line[valuestartpos] == '"' || line[valuestartpos] == '\'') {
+ if (line[valuestartpos] == '"' || line[valuestartpos] == '\'') {
quotechar = line[valuestartpos++];
// Scan until end of line or terminating quote
// valuestartpos++;
pos = valuestartpos;
- while(pos<line.length()) {
- if(line[pos++] == quotechar) { // Possibly end of quoted string
- if(line[pos] == quotechar) { // Protected quote
- pos++;
- } else {
- pos--;
- break;
- }
- }
+ while (pos < line.length()) {
+ if (line[pos++] == quotechar) { // Possibly end of quoted string
+ if (line[pos] == quotechar) { // Protected quote
+ pos++;
+ } else {
+ pos--;
+ break;
+ }
+ }
}
- if(TextList) {
- // Should check that a numerical assignment doesn't exist, but for
- // now, the same variable name can be used for strings and numbers
- string varnames(varname);
- AddString(varnames, line.substr(valuestartpos,pos-valuestartpos));
+ if (TextList) {
+ // Should check that a numerical assignment doesn't exist, but for
+ // now, the same variable name can be used for strings and numbers
+ string varnames(varname);
+ AddString(varnames, line.substr(valuestartpos, pos - valuestartpos));
}
continue;
}
- TString values((line.substr(valuestartpos)).c_str());
- TObjArray *vararr = values.Tokenize(",");
- Int_t nvals = vararr->GetLast()+1;
+ TString values((line.substr(valuestartpos)).c_str());
+ TObjArray* vararr = values.Tokenize(",");
+ Int_t nvals = vararr->GetLast() + 1;
- Int_t* ip=0;
- Double_t* fp=0;
+ Int_t* ip = 0;
+ Double_t* fp = 0;
// or expressions
- for(Int_t i=0;(ttype==0&&i<nvals);i++) {
- TString valstr = ((TObjString *)vararr->At(i))->GetString();
- if(valstr.IsFloat()) { // Is a number
- if(valstr.Contains(".") || valstr.Contains("e",TString::kIgnoreCase)) {
- ttype = 1;
- break;
- }
+ for (Int_t i = 0; (ttype == 0 && i < nvals); i++) {
+ TString valstr = ((TObjString*)vararr->At(i))->GetString();
+ if (valstr.IsFloat()) { // Is a number
+ if (valstr.Contains(".") || valstr.Contains("e", TString::kIgnoreCase)) {
+ ttype = 1;
+ break;
+ }
} else {
- ttype = 2; // Force float if expression or Var
- break;
+ ttype = 2; // Force float if expression or Var
+ break;
}
}
@@ -382,150 +379,145 @@ The ENGINE CTP support parameter "blocks" which were marked with
//
// There is some code duplication here. Refactor later
- Int_t newlength = currentindex + nvals;
- THaVar* existingvar=Find(varname);
- if(existingvar) {
+ Int_t newlength = currentindex + nvals;
+ THaVar* existingvar = Find(varname);
+ if (existingvar) {
string existingcomment;
existingcomment.assign(existingvar->GetTitle());
- if(!existingcomment.empty()) {
- current_comment.assign(existingcomment);
+ if (!existingcomment.empty()) {
+ current_comment.assign(existingcomment);
}
- Int_t existingtype=existingvar->GetType();
- Int_t existinglength=existingvar->GetLen();
- if(newlength > existinglength ||
- (existingtype == kInt && ttype > 0)) { // Length or type change needed
- if(newlength < existinglength) newlength = existinglength;
- Int_t newtype=-1;
- if(ttype>0 || existingtype == kDouble) {
- newtype = kDouble;
- fp = new Double_t[newlength];
- if(existingtype == kDouble) {
- Double_t* existingp= (Double_t*) existingvar->GetValuePointer();
- for(Int_t i=0;i<existinglength;i++) {
- fp[i] = existingp[i];
- }
- } else if(existingtype == kInt) {
- Int_t* existingp= (Int_t*) existingvar->GetValuePointer();
- for(Int_t i=0;i<existinglength;i++) {
- fp[i] = existingp[i];
- }
- } else {
- cout << "Whoops!" << endl;
- }
- } else if(existingtype == kInt) { // Stays int
- newtype = kInt;
- ip = new Int_t[newlength];
- Int_t* existingp= (Int_t*) existingvar->GetValuePointer();
- for(Int_t i=0;i<existinglength;i++) {
- ip[i] = existingp[i];
- }
- } else {
- cout << "Whoops!" << endl;
- }
- // Now copy new values in
- for(Int_t i=0;i<nvals;i++) {
- TString valstr = ((TObjString *)vararr->At(i))->GetString();
- if(newtype == kInt) {
- ip[currentindex+i] = valstr.Atoi();
- } else {
- if(valstr.IsFloat()) {
- fp[currentindex+i] = valstr.Atof();
- } else {
- THaFormula* formula = new THaFormula
- ("temp",valstr.Data(), (Bool_t) 0, this, 0);
- fp[currentindex+i] = formula->Eval();
- delete formula;
- }
- }
- }
- currentindex += nvals;
- // Remove old variable and recreate
- if(existingtype == kDouble) {
- delete [] (Double_t*) existingvar->GetValuePointer();
- } else if (existingtype == kInt) {
- delete [] (Int_t*) existingvar->GetValuePointer();
- }
- RemoveName(varname);
- char *arrayname=new char [strlen(varname)+20];
- sprintf(arrayname,"%s[%d]",varname,newlength);
- if(newtype == kInt) {
- Define(arrayname, current_comment.c_str(), *ip);
- } else {
- Define(arrayname, current_comment.c_str(), *fp);
- }
- delete[] arrayname;
+ Int_t existingtype = existingvar->GetType();
+ Int_t existinglength = existingvar->GetLen();
+ if (newlength > existinglength ||
+ (existingtype == kInt && ttype > 0)) { // Length or type change needed
+ if (newlength < existinglength)
+ newlength = existinglength;
+ Int_t newtype = -1;
+ if (ttype > 0 || existingtype == kDouble) {
+ newtype = kDouble;
+ fp = new Double_t[newlength];
+ if (existingtype == kDouble) {
+ Double_t* existingp = (Double_t*)existingvar->GetValuePointer();
+ for (Int_t i = 0; i < existinglength; i++) {
+ fp[i] = existingp[i];
+ }
+ } else if (existingtype == kInt) {
+ Int_t* existingp = (Int_t*)existingvar->GetValuePointer();
+ for (Int_t i = 0; i < existinglength; i++) {
+ fp[i] = existingp[i];
+ }
+ } else {
+ cout << "Whoops!" << endl;
+ }
+ } else if (existingtype == kInt) { // Stays int
+ newtype = kInt;
+ ip = new Int_t[newlength];
+ Int_t* existingp = (Int_t*)existingvar->GetValuePointer();
+ for (Int_t i = 0; i < existinglength; i++) {
+ ip[i] = existingp[i];
+ }
+ } else {
+ cout << "Whoops!" << endl;
+ }
+ // Now copy new values in
+ for (Int_t i = 0; i < nvals; i++) {
+ TString valstr = ((TObjString*)vararr->At(i))->GetString();
+ if (newtype == kInt) {
+ ip[currentindex + i] = valstr.Atoi();
+ } else {
+ if (valstr.IsFloat()) {
+ fp[currentindex + i] = valstr.Atof();
+ } else {
+ THaFormula* formula = new THaFormula("temp", valstr.Data(), (Bool_t)0, this, 0);
+ fp[currentindex + i] = formula->Eval();
+ delete formula;
+ }
+ }
+ }
+ currentindex += nvals;
+ // Remove old variable and recreate
+ if (existingtype == kDouble) {
+ delete[](Double_t*) existingvar->GetValuePointer();
+ } else if (existingtype == kInt) {
+ delete[](Int_t*) existingvar->GetValuePointer();
+ }
+ RemoveName(varname);
+ char* arrayname = new char[strlen(varname) + 20];
+ sprintf(arrayname, "%s[%d]", varname, newlength);
+ if (newtype == kInt) {
+ Define(arrayname, current_comment.c_str(), *ip);
+ } else {
+ Define(arrayname, current_comment.c_str(), *fp);
+ }
+ delete[] arrayname;
} else {
- // Existing array long enough and of right type, just copy to it.
- if(ttype == 0 && existingtype == kInt) {
- Int_t* existingp= (Int_t*) existingvar->GetValuePointer();
- for(Int_t i=0;i<nvals;i++) {
- TString valstr = ((TObjString *)vararr->At(i))->GetString();
- existingp[currentindex+i] = valstr.Atoi();
- }
- } else {
- Double_t* existingp= (Double_t*) existingvar->GetValuePointer();
- for(Int_t i=0;i<nvals;i++) {
- TString valstr = ((TObjString *)vararr->At(i))->GetString();
- if(valstr.IsFloat()) {
- existingp[currentindex+i] = valstr.Atof();
- } else {
- THaFormula* formula = new THaFormula
- ("temp",valstr.Data(), (Bool_t) 0, this, 0);
- existingp[currentindex+i] = formula->Eval();
- delete formula;
- }
- }
- }
- currentindex += nvals;
+ // Existing array long enough and of right type, just copy to it.
+ if (ttype == 0 && existingtype == kInt) {
+ Int_t* existingp = (Int_t*)existingvar->GetValuePointer();
+ for (Int_t i = 0; i < nvals; i++) {
+ TString valstr = ((TObjString*)vararr->At(i))->GetString();
+ existingp[currentindex + i] = valstr.Atoi();
+ }
+ } else {
+ Double_t* existingp = (Double_t*)existingvar->GetValuePointer();
+ for (Int_t i = 0; i < nvals; i++) {
+ TString valstr = ((TObjString*)vararr->At(i))->GetString();
+ if (valstr.IsFloat()) {
+ existingp[currentindex + i] = valstr.Atof();
+ } else {
+ THaFormula* formula = new THaFormula("temp", valstr.Data(), (Bool_t)0, this, 0);
+ existingp[currentindex + i] = formula->Eval();
+ delete formula;
+ }
+ }
+ }
+ currentindex += nvals;
}
} else {
- if(currentindex !=0) {
- cout << "currentindex=" << currentindex << " shouldn't be!" << endl;
+ if (currentindex != 0) {
+ cout << "currentindex=" << currentindex << " shouldn't be!" << endl;
}
- if(ttype==0) {
- ip = new Int_t[nvals];
+ if (ttype == 0) {
+ ip = new Int_t[nvals];
} else {
- fp = new Double_t[nvals];
+ fp = new Double_t[nvals];
}
- for(Int_t i=0;i<nvals;i++) {
- TString valstr = ((TObjString *)vararr->At(i))->GetString();
- if(ttype==0) {
- ip[i] = valstr.Atoi();
- } else {
- if(valstr.IsFloat()) {
- fp[i] = valstr.Atof();
- } else {
- THaFormula* formula = new THaFormula
- ("temp",valstr.Data(), (Bool_t) 0, this, 0);
- fp[i] = formula->Eval();
- delete formula;
- }
- }
+ for (Int_t i = 0; i < nvals; i++) {
+ TString valstr = ((TObjString*)vararr->At(i))->GetString();
+ if (ttype == 0) {
+ ip[i] = valstr.Atoi();
+ } else {
+ if (valstr.IsFloat()) {
+ fp[i] = valstr.Atof();
+ } else {
+ THaFormula* formula = new THaFormula("temp", valstr.Data(), (Bool_t)0, this, 0);
+ fp[i] = formula->Eval();
+ delete formula;
+ }
+ }
}
currentindex = nvals;
- char *arrayname=new char [strlen(varname)+20];
- sprintf(arrayname,"%s[%d]",varname,nvals);
- if(ttype==0) {
- Define(arrayname, current_comment.c_str(), *ip);
+ char* arrayname = new char[strlen(varname) + 20];
+ sprintf(arrayname, "%s[%d]", varname, nvals);
+ if (ttype == 0) {
+ Define(arrayname, current_comment.c_str(), *ip);
} else {
- Define(arrayname, current_comment.c_str(), *fp);
+ Define(arrayname, current_comment.c_str(), *fp);
}
delete[] arrayname;
}
- delete vararr; // Discard result of Tokenize
+ delete vararr; // Discard result of Tokenize
// cout << line << endl;
-
}
return;
-
}
//_____________________________________________________________________________
-Int_t THcParmList::LoadParmValues(const DBRequest* list, const char* prefix)
-{
+Int_t THcParmList::LoadParmValues(const DBRequest* list, const char* prefix) {
/**
The following example loads several parameters held in the `gHcParms`
parameter cache into scalar variables or arrays.
@@ -550,77 +542,80 @@ zero), then there will be no error if the parameter is missing.
*/
- const DBRequest *ti = list;
- Int_t cnt=0;
- Int_t this_cnt=0;
+ const DBRequest* ti = list;
+ Int_t cnt = 0;
+ Int_t this_cnt = 0;
- if( !prefix ) prefix = "";
+ if (!prefix)
+ prefix = "";
- while ( ti && ti->name ) {
- string keystr(prefix); keystr.append(ti->name);
+ while (ti && ti->name) {
+ string keystr(prefix);
+ keystr.append(ti->name);
const char* key = keystr.c_str();
// cout <<"Now at "<<ti->name<<endl;
this_cnt = 0;
- if(this->Find(key)) {
+ if (this->Find(key)) {
VarType ty = this->Find(key)->GetType();
- if (ti->nelem>1) {
- // it is an array, use the appropriateinterface
- switch (ti->type) {
- case (kDouble) :
- this_cnt = GetArray(key,static_cast<Double_t*>(ti->var),ti->nelem);
- break;
- case (kInt) :
- this_cnt = GetArray(key,static_cast<Int_t*>(ti->var),ti->nelem);
- break;
- default:
- Error("THcParmList","Invalid type to read %s",key);
- break;
- }
+ if (ti->nelem > 1) {
+ // it is an array, use the appropriateinterface
+ switch (ti->type) {
+ case (kDouble):
+ this_cnt = GetArray(key, static_cast<Double_t*>(ti->var), ti->nelem);
+ break;
+ case (kInt):
+ this_cnt = GetArray(key, static_cast<Int_t*>(ti->var), ti->nelem);
+ break;
+ default:
+ Error("THcParmList", "Invalid type to read %s", key);
+ break;
+ }
} else {
- switch (ti->type) {
- case (kDouble) :
- if(ty == kInt) {
- *static_cast<Double_t*>(ti->var)=*(Int_t *)this->Find(key)->GetValuePointer();
- } else if (ty == kDouble) {
- *static_cast<Double_t*>(ti->var)=*(Double_t *)this->Find(key)->GetValuePointer();
- } else {
- cout << "*** ERROR!!! Type Mismatch " << key << endl;
- }
- this_cnt=1;
-
- break;
- case (kInt) :
- if(ty == kInt) {
- *static_cast<Int_t*>(ti->var)=*(Int_t *)this->Find(key)->GetValuePointer();
- } else if (ty == kDouble) {
- *static_cast<Int_t*>(ti->var)=TMath::Nint(*(Double_t *)this->Find(key)->GetValuePointer());
- _logger->warn("Rounded {} to nearest integer ",key);
- } else {
- _logger->error("Type Mismatch {} ",key);
- }
- this_cnt=1;
- break;
- default:
- _logger->error("THcParmList: Invalid type to read {} ",key);
+ switch (ti->type) {
+ case (kDouble):
+ if (ty == kInt) {
+ *static_cast<Double_t*>(ti->var) = *(Int_t*)this->Find(key)->GetValuePointer();
+ } else if (ty == kDouble) {
+ *static_cast<Double_t*>(ti->var) = *(Double_t*)this->Find(key)->GetValuePointer();
+ } else {
+ cout << "*** ERROR!!! Type Mismatch " << key << endl;
+ }
+ this_cnt = 1;
+
break;
- }
+ case (kInt):
+ if (ty == kInt) {
+ *static_cast<Int_t*>(ti->var) = *(Int_t*)this->Find(key)->GetValuePointer();
+ } else if (ty == kDouble) {
+ *static_cast<Int_t*>(ti->var) =
+ TMath::Nint(*(Double_t*)this->Find(key)->GetValuePointer());
+ _logger->warn("Rounded {} to nearest integer ", key);
+ } else {
+ _logger->error("Type Mismatch {} ", key);
+ }
+ this_cnt = 1;
+ break;
+ default:
+ _logger->error("THcParmList: Invalid type to read {} ", key);
+ break;
+ }
}
- } else { // See if it is a text variable
+ } else { // See if it is a text variable
const char* value = GetString(key);
- if(value) {
- this_cnt = 1;
- if(ti->type == kString) {
- *((string*)ti->var) = string(value);
- } else if (ti->type == kTString) {
- *((TString*)ti->var) = (TString) value;
- } else {
- Error("THcParmList","No conversion for strings: %s",key);
- }
+ if (value) {
+ this_cnt = 1;
+ if (ti->type == kString) {
+ *((string*)ti->var) = string(value);
+ } else if (ti->type == kTString) {
+ *((TString*)ti->var) = (TString)value;
+ } else {
+ Error("THcParmList", "No conversion for strings: %s", key);
+ }
}
}
- if (this_cnt<=0) {
- if ( !ti->optional ) {
+ if (this_cnt <= 0) {
+ if (!ti->optional) {
string msg = string("Could not find `") + key + "` in database!";
throw std::runtime_error("<THcParmList::LoadParmValues>: " + msg);
}
@@ -633,67 +628,68 @@ zero), then there will be no error if the parameter is missing.
// READING AN ARRAY INTO A C-style ARRAY
//_____________________________________________________________________________
-Int_t THcParmList::GetArray(const char* attr, Int_t* array, Int_t size)
-{
- return ReadArray(attr,array,size);
+Int_t THcParmList::GetArray(const char* attr, Int_t* array, Int_t size) {
+ return ReadArray(attr, array, size);
}
//_____________________________________________________________________________
-Int_t THcParmList::GetArray(const char* attr, Double_t* array, Int_t size)
-{
- return ReadArray(attr,array,size);
+Int_t THcParmList::GetArray(const char* attr, Double_t* array, Int_t size) {
+ return ReadArray(attr, array, size);
}
//_____________________________________________________________________________
-template<class T>
-Int_t THcParmList::ReadArray(const char* attrC, T* array, Int_t size)
-{
+template <class T>
+Int_t THcParmList::ReadArray(const char* attrC, T* array, Int_t size) {
/**
No resizing is done, so only 'size' elements may be stored.
*/
- Int_t cnt=0;
+ Int_t cnt = 0;
- THaVar *var = Find(attrC);
- if(!var) return(cnt);
+ THaVar* var = Find(attrC);
+ if (!var)
+ return (cnt);
VarType ty = var->GetType();
- if( ty != kInt && ty != kDouble) {
+ if (ty != kInt && ty != kDouble) {
cout << "*** ERROR: " << attrC << " has unsupported type " << ty << endl;
- return(cnt);
+ return (cnt);
}
- Int_t sz = var->GetLen();
- const void *vp = var->GetValuePointer();
- if(size != sz) {
+ Int_t sz = var->GetLen();
+ const void* vp = var->GetValuePointer();
+ if (size > sz) {
+ _logger->error("requested {} elements of {} which has only {} elements", size, attrC, sz);
+ exit(EXIT_FAILURE);
+ } else if (size < sz) {
_logger->warn("requested {} elements of {} which has length {}", size, attrC, sz);
}
- if(size<sz) sz = size;
+ if (size < sz)
+ sz = size;
Int_t donint = 0;
- if(ty == kDouble && typeid(array[0]) == typeid(Int_t)) {
- donint = 1; // Use nint when putting doubles in nint
+ if (ty == kDouble && typeid(array[0]) == typeid(Int_t)) {
+ donint = 1; // Use nint when putting doubles in nint
cout << "*** WARNING!!! Rounded " << attrC << " elements to nearest integer " << endl;
}
- for(cnt=0;cnt<sz;cnt++) {
- if(ty == kInt) {
+ for (cnt = 0; cnt < sz; cnt++) {
+ if (ty == kInt) {
array[cnt] = ((Int_t*)vp)[cnt];
- } else
- if(donint) {
- array[cnt] = TMath::Nint(((Double_t*)vp)[cnt]);
- } else {
- array[cnt] = ((Double_t*)vp)[cnt];
- }
+ } else if (donint) {
+ array[cnt] = TMath::Nint(((Double_t*)vp)[cnt]);
+ } else {
+ array[cnt] = ((Double_t*)vp)[cnt];
+ }
}
- return(cnt);
+ return (cnt);
}
//_____________________________________________________________________________
-std::string THcParmList::PrintJSON(int run_number ) const {
- TIter next(this);
+std::string THcParmList::PrintJSON(int run_number) const {
+ TIter next(this);
nlohmann::json j;
- while( THaVar* obj = (THaVar*) next() ) {
+ while (THaVar* obj = (THaVar*)next()) {
- if( obj->IsBasic() ) {
- if( obj->IsArray() ) {
- if( obj->GetLen() == 1 ) {
+ if (obj->IsBasic()) {
+ if (obj->IsArray()) {
+ if (obj->GetLen() == 1) {
j[obj->GetName()] = obj->GetValue();
} else {
j[obj->GetName()] = obj->GetValues();
@@ -701,133 +697,127 @@ std::string THcParmList::PrintJSON(int run_number ) const {
} else {
obj->Print();
}
- } else {
+ } else {
obj->Print();
}
}
nlohmann::json jrun;
jrun[std::to_string(run_number)] = j;
- //std::cout << j.dump(2) << "\n";
+ // std::cout << j.dump(2) << "\n";
// write prettified JSON to another file
- //std::ofstream o("pretty.json");
- //o << std::setw(4) << jrun << std::endl;
+ // std::ofstream o("pretty.json");
+ // o << std::setw(4) << jrun << std::endl;
return jrun.dump(2);
}
-
-void THcParmList::PrintFull( Option_t* option ) const
-{
+void THcParmList::PrintFull(Option_t* option) const {
THaVarList::PrintFull(option);
TextList->Print();
}
#ifdef WITH_CCDB
//_____________________________________________________________________________
-Int_t THcParmList::OpenCCDB(Int_t runnum)
-{
+Int_t THcParmList::OpenCCDB(Int_t runnum) {
// Use the Environment variable "CCDB_CONNECTION" as the
// connection string
const char* connection_string = gSystem->Getenv("CCDB_CONNECTION");
- return(OpenCCDB(runnum,connection_string));
+ return (OpenCCDB(runnum, connection_string));
}
-Int_t THcParmList::OpenCCDB(Int_t runnum, const char* connection_string)
-{
+Int_t THcParmList::OpenCCDB(Int_t runnum, const char* connection_string) {
// Connect to a CCDB database pointed to by connection_string
- std::string s (connection_string);
- CCDB_obj = new SQLiteCalibration(runnum);
+ std::string s(connection_string);
+ CCDB_obj = new SQLiteCalibration(runnum);
Int_t result = CCDB_obj->Connect(s);
- if(!result) return -1; // Need some error codes
+ if (!result)
+ return -1; // Need some error codes
cout << "Opened " << s << " for run " << runnum << endl;
return 0;
}
-Int_t THcParmList::CloseCCDB()
-{
+Int_t THcParmList::CloseCCDB() {
delete CCDB_obj;
- return(0);
+ return (0);
}
-Int_t THcParmList::LoadCCDBDirectory(const char* directory,
- const char* prefix)
-{
+Int_t THcParmList::LoadCCDBDirectory(const char* directory, const char* prefix) {
// Load all parameters in directory
// Prepend prefix onto the name of each
- std::string dirname (directory);
+ std::string dirname(directory);
- if(dirname[dirname.length()-1]!='/') {
+ if (dirname[dirname.length() - 1] != '/') {
dirname.append("/");
}
- Int_t dirlen=dirname.length();
+ Int_t dirlen = dirname.length();
vector<string> namepaths;
CCDB_obj->GetListOfNamepaths(namepaths);
- for(UInt_t iname=0;iname<namepaths.size();iname++) {
- std::string varname (namepaths[iname]);
- if(varname.compare(0,dirlen,dirname) == 0) {
- varname.replace(0,dirlen,prefix);
+ for (UInt_t iname = 0; iname < namepaths.size(); iname++) {
+ std::string varname(namepaths[iname]);
+ if (varname.compare(0, dirlen, dirname) == 0) {
+ varname.replace(0, dirlen, prefix);
// cout << namepaths[iname] << " -> " << varname << endl;
// To what extent is there duplication here with Load() method?
// Retrieve assignment
- Assignment* assignment = CCDB_obj->GetAssignment(namepaths[iname], true);
- ConstantsTypeColumn::ColumnTypes ccdbtype=assignment->GetValueType(0);
- Int_t ccdbncolumns=assignment->GetColumnsCount();
- Int_t ccdbnrows=assignment->GetRowsCount();
- std::string title = assignment->GetTypeTable()->GetComment();
+ Assignment* assignment = CCDB_obj->GetAssignment(namepaths[iname], true);
+ ConstantsTypeColumn::ColumnTypes ccdbtype = assignment->GetValueType(0);
+ Int_t ccdbncolumns = assignment->GetColumnsCount();
+ Int_t ccdbnrows = assignment->GetRowsCount();
+ std::string title = assignment->GetTypeTable()->GetComment();
// Only load single column tables
- if(ccdbncolumns == 1) {
-
- THaVar* existingvar=Find(varname.c_str());
- // Need to append [size] to end of varname
- char sizestring[20];
- sprintf(sizestring,"[%d]",ccdbnrows);
- std::string size_str (sizestring);
- std::string varnamearray (varname);
- varnamearray.append(size_str);
-
- // Select data type
- if(ccdbtype==ConstantsTypeColumn::cIntColumn) {
- vector<vector<int> > data;
- CCDB_obj->GetCalib(data, namepaths[iname]);
-
- if(existingvar) {
- RemoveName(varname.c_str());
- }
-
- Int_t* ip = new Int_t[data.size()];
- for(UInt_t row=0;row<data.size(); row++) {
- ip[row] = data[row][0];
- }
- Define(varnamearray.c_str(), title.c_str(), *ip);
-
- } else if (ccdbtype==ConstantsTypeColumn::cDoubleColumn) {
- vector<vector<double> > data;
- CCDB_obj->GetCalib(data, namepaths[iname]);
-
- if(existingvar) {
- RemoveName(varname.c_str());
- }
-
- Double_t* fp = new Double_t[data.size()];
- for(UInt_t row=0;row<data.size(); row++) {
- fp[row] = data[row][0];
- }
- Define(varnamearray.c_str(), title.c_str(), *fp);
- } else if (ccdbtype==ConstantsTypeColumn::cStringColumn) {
- if(ccdbnrows > 1) {
- cout << namepaths[iname] << ": Only first element of CCDB string array loaded." << endl;
- }
- vector<vector<string> > data;
- CCDB_obj->GetCalib(data, namepaths[iname]);
- AddString(varname, data[0][0]);
- } else {
- cout << namepaths[iname] << ": Unsupported CCDB data type: " << ccdbtype << endl;
- }
+ if (ccdbncolumns == 1) {
+
+ THaVar* existingvar = Find(varname.c_str());
+ // Need to append [size] to end of varname
+ char sizestring[20];
+ sprintf(sizestring, "[%d]", ccdbnrows);
+ std::string size_str(sizestring);
+ std::string varnamearray(varname);
+ varnamearray.append(size_str);
+
+ // Select data type
+ if (ccdbtype == ConstantsTypeColumn::cIntColumn) {
+ vector<vector<int>> data;
+ CCDB_obj->GetCalib(data, namepaths[iname]);
+
+ if (existingvar) {
+ RemoveName(varname.c_str());
+ }
+
+ Int_t* ip = new Int_t[data.size()];
+ for (UInt_t row = 0; row < data.size(); row++) {
+ ip[row] = data[row][0];
+ }
+ Define(varnamearray.c_str(), title.c_str(), *ip);
+
+ } else if (ccdbtype == ConstantsTypeColumn::cDoubleColumn) {
+ vector<vector<double>> data;
+ CCDB_obj->GetCalib(data, namepaths[iname]);
+
+ if (existingvar) {
+ RemoveName(varname.c_str());
+ }
+
+ Double_t* fp = new Double_t[data.size()];
+ for (UInt_t row = 0; row < data.size(); row++) {
+ fp[row] = data[row][0];
+ }
+ Define(varnamearray.c_str(), title.c_str(), *fp);
+ } else if (ccdbtype == ConstantsTypeColumn::cStringColumn) {
+ if (ccdbnrows > 1) {
+ cout << namepaths[iname] << ": Only first element of CCDB string array loaded." << endl;
+ }
+ vector<vector<string>> data;
+ CCDB_obj->GetCalib(data, namepaths[iname]);
+ AddString(varname, data[0][0]);
+ } else {
+ cout << namepaths[iname] << ": Unsupported CCDB data type: " << ccdbtype << endl;
+ }
} else {
- cout << namepaths[iname] << ": Multicolumn CCDB variables not supported" << endl;
+ cout << namepaths[iname] << ": Multicolumn CCDB variables not supported" << endl;
}
}
}
diff --git a/src/THcScalerEvtHandler.cxx b/src/THcScalerEvtHandler.cxx
index 1f06497096487779ffd852c868f87086d6b2b366..19a4a6dd4ee2f31d77aed31a381690e428caa499 100644
--- a/src/THcScalerEvtHandler.cxx
+++ b/src/THcScalerEvtHandler.cxx
@@ -40,7 +40,6 @@ To enable debugging you may try this in the setup script
#include <map>
#include <sstream>
-#include "THcScalerEvtHandler.h"
#include "GenScaler.h"
#include "Helper.h"
#include "Scaler1151.h"
@@ -56,6 +55,7 @@ To enable debugging you may try this in the setup script
#include "THaVarList.h"
#include "THcGlobals.h"
#include "THcParmList.h"
+#include "THcScalerEvtHandler.h"
#include "TMath.h"
#include "TNamed.h"
#include "TROOT.h"
@@ -151,7 +151,7 @@ Int_t THcScalerEvtHandler::ReadDatabase(const TDatime& date) {
fBCM_SatQuadratic = new Double_t[fNumBCMs];
fBCM_delta_charge = new Double_t[fNumBCMs];
string bcm_namelist;
- DBRequest list2[] = {{"BCM_Gain", fBCM_Gain, kDouble, (UInt_t)fNumBCMs},
+ DBRequest list2[] = {{"BCM_Gain", fBCM_Gain, kDouble, (UInt_t)fNumBCMs},
{"BCM_Offset", fBCM_Offset, kDouble, (UInt_t)fNumBCMs},
{"BCM_SatQuadratic", fBCM_SatQuadratic, kDouble, (UInt_t)fNumBCMs, 1},
{"BCM_SatOffset", fBCM_SatOffset, kDouble, (UInt_t)fNumBCMs, 1},
@@ -383,6 +383,7 @@ Int_t THcScalerEvtHandler::AnalyzeBuffer(UInt_t* rdata, Bool_t onlysync) {
size_t ivar = scalerloc[i]->ivar;
size_t idx = scalerloc[i]->index;
size_t ichan = scalerloc[i]->ichan;
+
if (evcount == 0) {
if (fDebugFile)
*fDebugFile << "Debug dvarsFirst " << i << " " << ivar << " " << idx << " " << ichan
@@ -421,7 +422,8 @@ Int_t THcScalerEvtHandler::AnalyzeBuffer(UInt_t* rdata, Bool_t onlysync) {
fBCM_Gain[bcm_ind];
dvars[ivar] =
dvars[ivar] +
- fBCM_SatOffset[bcm_ind] * TMath::Max(dvars[ivar] - fBCM_SatOffset[i], 0.0);
+ fBCM_SatQuadratic[bcm_ind] *
+ TMath::Power(TMath::Max(dvars[ivar] - fBCM_SatOffset[bcm_ind], 0.0), 2.0);
}
if (bcm_ind == fbcm_Current_Threshold_Index)
scal_current = dvars[ivar];
@@ -478,7 +480,7 @@ Int_t THcScalerEvtHandler::AnalyzeBuffer(UInt_t* rdata, Bool_t onlysync) {
dvarsFirst[ivar] =
dvarsFirst[ivar] +
fBCM_SatQuadratic[bcm_ind] *
- TMath::Power(TMath::Max(dvars[ivar] - fBCM_SatOffset[i], 0.0), 2.);
+ TMath::Power(TMath::Max(dvars[ivar] - fBCM_SatOffset[bcm_ind], 0.0), 2.);
}
if (bcm_ind == fbcm_Current_Threshold_Index)
scal_current = dvarsFirst[ivar];
@@ -504,6 +506,7 @@ Int_t THcScalerEvtHandler::AnalyzeBuffer(UInt_t* rdata, Bool_t onlysync) {
} else {
cout << "THcScalerEvtHandler:: ERROR:: incorrect index " << ivar << " " << idx << " "
<< ichan << endl;
+>>>>>>> c420edd53c56a5a9d8db08d317c7f15c82de6d17
}
} else { // evcount != 0
if (fDebugFile)
@@ -588,7 +591,7 @@ Int_t THcScalerEvtHandler::AnalyzeBuffer(UInt_t* rdata, Bool_t onlysync) {
}
// printf("event %i index %i fBCMname %s scalerloc %s offset %f gain %f computed
//%f\n",evcount, bcm_ind,
- //fBCM_Name[bcm_ind],scalerloc[ivar]->name.Data(),fBCM_Offset[bcm_ind],fBCM_Gain[bcm_ind],dvars[ivar]);
+ // fBCM_Name[bcm_ind],scalerloc[ivar]->name.Data(),fBCM_Offset[bcm_ind],fBCM_Gain[bcm_ind],dvars[ivar]);
}
if (fDebugFile)
*fDebugFile << " dvars " << scalerloc[ivar]->ikind << " " << dvars[ivar] << endl;
@@ -638,7 +641,7 @@ Int_t THcScalerEvtHandler::AnalyzeBuffer(UInt_t* rdata, Bool_t onlysync) {
evcount = evcount + 1;
evcountR = evcount;
//
- for (size_t j = 0; j < scal_prev_read.size(); j++)
+ for (size_t j = 0; j < scal_prev_read.size(); j++)
scal_prev_read[j] = scal_present_read[j];
//
for (size_t j = 0; j < scalers.size(); j++)
@@ -703,7 +706,7 @@ THaAnalysisObject::EStatus THcScalerEvtHandler::Init(const TDatime& date) {
if (pos1 != minus1 && dbline.size() > 4) {
string sdesc = "";
for (size_t j = 5; j < dbline.size(); j++)
- sdesc = sdesc + " " + dbline[j];
+ sdesc = sdesc + " " + dbline[j];
UInt_t islot = atoi(dbline[1].c_str());
UInt_t ichan = atoi(dbline[2].c_str());
UInt_t ikind = atoi(dbline[3].c_str());
diff --git a/src/THcShowerHit.h b/src/THcShowerHit.h
index 7619fb081b27e7a04e204f71cce4ae274977c992..41bfd17c24d4068d3ad8fd8a5961feb858ba9d03 100644
--- a/src/THcShowerHit.h
+++ b/src/THcShowerHit.h
@@ -7,6 +7,7 @@
#include <iterator>
#include <iostream>
#include <memory>
+#include <vector>
#include "TMath.h"
using namespace std;
diff --git a/tests/README.md b/tests/README.md
deleted file mode 100644
index 6a0afb3ebdafedc26379432edcae7ee373367894..0000000000000000000000000000000000000000
--- a/tests/README.md
+++ /dev/null
@@ -1,13 +0,0 @@
-Hcana tests
-===========
-
-## Current tests:
-
-- ep elastic tests
-
-## Tests to add:
-
-- Start time check
-- Focal plane times test
-- "Tracking efficiency" test
-- Invariant mass check (jpsi)
diff --git a/tests/elastic_coin_replay.cxx b/tests/elastic_coin_replay.cxx
deleted file mode 100644
index 7011fa8a309be4d7572aaf931f7ac7df1496ea50..0000000000000000000000000000000000000000
--- a/tests/elastic_coin_replay.cxx
+++ /dev/null
@@ -1,296 +0,0 @@
-#include <fmt/core.h>
-#include <fmt/ostream.h>
-#include <vector>
-
-#include "TString.h"
-
-R__LOAD_LIBRARY(libHallC.so)
-#include "hcana/HallC_Data.h"
-#include "DecData.h"
-//R__LOAD_LIBRARY(libScandalizer.so)
-//#include "monitor/DetectorDisplay.h"
-//#include "monitor/DisplayPlots.h"
-//#include "monitor/MonitoringDisplay.h"
-//#include "scandalizer/PostProcessors.h"
-//#include "scandalizer/ScriptHelpers.h"
-//
-//#include "THaPostProcess.h"
-//#include "monitor/ExperimentMonitor.h"
-//#include "scandalizer/PostProcessors.h"
-#include "THcAnalyzer.h"
-#include "THaCut.h"
-#include "THcGlobals.h"
-#include "THcHallCSpectrometer.h"
-#include "THcDetectorMap.h"
-#include "THcCherenkov.h"
-#include "THcDC.h"
-#include "THcHodoscope.h"
-#include "THcParmList.h"
-#include "THaGoldenTrack.h"
-#include "THcHodoEff.h"
-#include "THcScalerEvtHandler.h"
-#include "THcShower.h"
-#include "THcReactionPoint.h"
-#include "THcExtTarCor.h"
-#include "THcRasteredBeam.h"
-#include "THcRun.h"
-#include "THcCoinTime.h"
-#include "THcConfigEvtHandler.h"
-#include "THcTrigDet.h"
-#include "THcTrigApp.h"
-#include "THcSecondaryKine.h"
-#include "THcAerogel.h"
-#include "THcPrimaryKine.h"
-#include "THaReactionPoint.h"
-
-
-
-void elastic_coin_replay(Int_t RunNumber = 0, Int_t MaxEvent = -1) {
- using namespace std;
-
- // Get RunNumber and MaxEvent if not provided.
- if( RunNumber<=0 ) {
- std::exit(-1);
- }
-
- // Create file name patterns.
- const char* RunFileNamePattern = "coin_all_%05d.dat";
- vector<TString> pathList;
- pathList.push_back(".");
- pathList.push_back("./raw");
- pathList.push_back("./raw/../raw.copiedtotape");
- pathList.push_back("./cache");
-
- //const char* RunFileNamePattern = "raw/coin_all_%05d.dat";
- const char* ROOTFileNamePattern = "ROOTfiles/coin_replay_production_%d_%d.root";
-
- // Load global parameters
- gHcParms->Define("gen_run_number", "Run Number", RunNumber);
- gHcParms->AddString("g_ctp_database_filename", "DBASE/COIN/standard.database");
- gHcParms->Load(gHcParms->GetString("g_ctp_database_filename"), RunNumber);
- gHcParms->Load(gHcParms->GetString("g_ctp_parm_filename"));
- gHcParms->Load(gHcParms->GetString("g_ctp_kinematics_filename"), RunNumber);
- // Load params for COIN trigger configuration
- gHcParms->Load("PARAM/TRIG/tcoin.param");
- // Load fadc debug parameters
- gHcParms->Load("PARAM/HMS/GEN/h_fadc_debug.param");
- gHcParms->Load("PARAM/SHMS/GEN/p_fadc_debug.param");
-
- // const char* CurrentFileNamePattern = "low_curr_bcm/bcmcurrent_%d.param";
- // gHcParms->Load(Form(CurrentFileNamePattern, RunNumber));
-
- // Load the Hall C detector map
- gHcDetectorMap = new THcDetectorMap();
- gHcDetectorMap->Load("MAPS/COIN/DETEC/coin.map");
-
- // Dec data
- gHaApps->Add(new Podd::DecData("D","Decoder raw data"));
- //=:=:=:=
- // SHMS
- //=:=:=:=
-
- // Set up the equipment to be analyzed.
- THcHallCSpectrometer* SHMS = new THcHallCSpectrometer("P", "SHMS");
- SHMS->SetEvtType(1);
- SHMS->AddEvtType(4);
- SHMS->AddEvtType(5);
- SHMS->AddEvtType(6);
- SHMS->AddEvtType(7);
- gHaApps->Add(SHMS);
- // Add Noble Gas Cherenkov to SHMS apparatus
- THcCherenkov* pngcer = new THcCherenkov("ngcer", "Noble Gas Cherenkov");
- SHMS->AddDetector(pngcer);
- // Add drift chambers to SHMS apparatus
- THcDC* pdc = new THcDC("dc", "Drift Chambers");
- SHMS->AddDetector(pdc);
- // Add hodoscope to SHMS apparatus
- THcHodoscope* phod = new THcHodoscope("hod", "Hodoscope");
- SHMS->AddDetector(phod);
- // Add Heavy Gas Cherenkov to SHMS apparatus
- THcCherenkov* phgcer = new THcCherenkov("hgcer", "Heavy Gas Cherenkov");
- SHMS->AddDetector(phgcer);
- // Add Aerogel Cherenkov to SHMS apparatus
- THcAerogel* paero = new THcAerogel("aero", "Aerogel");
- SHMS->AddDetector(paero);
- // Add calorimeter to SHMS apparatus
- THcShower* pcal = new THcShower("cal", "Calorimeter");
- SHMS->AddDetector(pcal);
-
- // THcBCMCurrent* hbc = new THcBCMCurrent("H.bcm", "BCM current check");
- // gHaPhysics->Add(hbc);
-
- // Add rastered beam apparatus
- THaApparatus* pbeam = new THcRasteredBeam("P.rb", "Rastered Beamline");
- gHaApps->Add(pbeam);
- // Add physics modules
- // Calculate reaction point
- THcReactionPoint* prp = new THcReactionPoint("P.react", "SHMS reaction point", "P", "P.rb");
- gHaPhysics->Add(prp);
- // Calculate extended target corrections
- THcExtTarCor* pext = new THcExtTarCor("P.extcor", "HMS extended target corrections", "P", "P.react");
- gHaPhysics->Add(pext);
- // Calculate golden track quantites
- THaGoldenTrack* pgtr = new THaGoldenTrack("P.gtr", "SHMS Golden Track", "P");
- gHaPhysics->Add(pgtr);
- // Calculate the hodoscope efficiencies
- THcHodoEff* peff = new THcHodoEff("phodeff", "SHMS hodo efficiency", "P.hod");
- gHaPhysics->Add(peff);
-
- // Add event handler for scaler events
- THcScalerEvtHandler* pscaler = new THcScalerEvtHandler("P", "Hall C scaler event type 1");
- pscaler->AddEvtType(1);
- pscaler->AddEvtType(4);
- pscaler->AddEvtType(5);
- pscaler->AddEvtType(6);
- pscaler->AddEvtType(7);
- pscaler->AddEvtType(129);
- pscaler->SetDelayedType(129);
- pscaler->SetUseFirstEvent(kTRUE);
- gHaEvtHandlers->Add(pscaler);
-
- //=:=:=
- // HMS
- //=:=:=
-
- // Set up the equipment to be analyzed.
- THcHallCSpectrometer* HMS = new THcHallCSpectrometer("H", "HMS");
- HMS->SetEvtType(2);
- HMS->AddEvtType(4);
- HMS->AddEvtType(5);
- HMS->AddEvtType(6);
- HMS->AddEvtType(7);
- gHaApps->Add(HMS);
- // Add drift chambers to HMS apparatus
- THcDC* hdc = new THcDC("dc", "Drift Chambers");
- HMS->AddDetector(hdc);
- // Add hodoscope to HMS apparatus
- THcHodoscope* hhod = new THcHodoscope("hod", "Hodoscope");
- HMS->AddDetector(hhod);
- // Add Cherenkov to HMS apparatus
- THcCherenkov* hcer = new THcCherenkov("cer", "Heavy Gas Cherenkov");
- HMS->AddDetector(hcer);
- // Add Aerogel Cherenkov to HMS apparatus
- // THcAerogel* haero = new THcAerogel("aero", "Aerogel");
- // HMS->AddDetector(haero);
- // Add calorimeter to HMS apparatus
- THcShower* hcal = new THcShower("cal", "Calorimeter");
- HMS->AddDetector(hcal);
-
- // Add rastered beam apparatus
- THaApparatus* hbeam = new THcRasteredBeam("H.rb", "Rastered Beamline");
- gHaApps->Add(hbeam);
- // Add physics modules
- // Calculate reaction point
- THcReactionPoint* hrp = new THcReactionPoint("H.react", "HMS reaction point", "H", "H.rb");
- gHaPhysics->Add(hrp);
- // Calculate extended target corrections
- THcExtTarCor* hext = new THcExtTarCor("H.extcor", "HMS extended target corrections", "H", "H.react");
- gHaPhysics->Add(hext);
- // Calculate golden track quantities
- THaGoldenTrack* hgtr = new THaGoldenTrack("H.gtr", "HMS Golden Track", "H");
- gHaPhysics->Add(hgtr);
- // Calculate the hodoscope efficiencies
- THcHodoEff* heff = new THcHodoEff("hhodeff", "HMS hodo efficiency", "H.hod");
- gHaPhysics->Add(heff);
-
- // Add event handler for scaler events
- THcScalerEvtHandler *hscaler = new THcScalerEvtHandler("H", "Hall C scaler event type 4");
- hscaler->AddEvtType(2);
- hscaler->AddEvtType(4);
- hscaler->AddEvtType(5);
- hscaler->AddEvtType(6);
- hscaler->AddEvtType(7);
- hscaler->AddEvtType(129);
- hscaler->SetDelayedType(129);
- hscaler->SetUseFirstEvent(kTRUE);
- gHaEvtHandlers->Add(hscaler);
-
- //=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=
- // Kinematics Modules
- //=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=
-
- // Add Physics Module to calculate primary (scattered electrons) beam kinematics
- THcPrimaryKine* hkin_primary = new THcPrimaryKine("H.kin.primary", "HMS Single Arm Kinematics", "H", "H.rb");
- gHaPhysics->Add(hkin_primary);
- // Add Physics Module to calculate secondary (scattered hadrons) beam kinematics
- THcSecondaryKine* pkin_secondary = new THcSecondaryKine("P.kin.secondary", "SHMS Single Arm Kinematics", "P", "H.kin.primary");
- gHaPhysics->Add(pkin_secondary);
-
- //=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=
- // Global Objects & Event Handlers
- //=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=:=
-
- // Add trigger apparatus
- THaApparatus* TRG = new THcTrigApp("T", "TRG");
- gHaApps->Add(TRG);
- // Add trigger detector to trigger apparatus
- THcTrigDet* coin = new THcTrigDet("coin", "Coincidence Trigger Information");
- // Suppress missing reference time warnings for these event types
- coin->SetEvtType(1);
- coin->AddEvtType(2);
- TRG->AddDetector(coin);
-
-
- //Add coin physics module THcCoinTime::THcCoinTime (const char *name, const char* description, const char* hadArmName,
- // const char* elecArmName, const char* coinname) :
- THcCoinTime* coinTime = new THcCoinTime("CTime", "Coincidende Time Determination", "P", "H", "T.coin");
- gHaPhysics->Add(coinTime);
-
- // Add event handler for prestart event 125.
- THcConfigEvtHandler* ev125 = new THcConfigEvtHandler("HC", "Config Event type 125");
- gHaEvtHandlers->Add(ev125);
- // Add event handler for EPICS events
- THaEpicsEvtHandler* hcepics = new THaEpicsEvtHandler("epics", "HC EPICS event type 180");
- gHaEvtHandlers->Add(hcepics);
-
- // Set up the analyzer - we use the standard one,
- // but this could be an experiment-specific one as well.
- // The Analyzer controls the reading of the data, executes
- // tests/cuts, loops over Acpparatus's and PhysicsModules,
- // and executes the output routines.
- THcAnalyzer* analyzer = new THcAnalyzer;
-
- // A simple event class to be output to the resulting tree.
- // Creating your own descendant of THaEvent is one way of
- // defining and controlling the output.
- THaEvent* event = new THaEvent;
-
- // Define the run(s) that we want to analyze.
- // We just set up one, but this could be many.
- THcRun* run = new THcRun( pathList, Form(RunFileNamePattern, RunNumber) );
-
- // Set to read in Hall C run database parameters
- run->SetRunParamClass("THcRunParameters");
-
- // Eventually need to learn to skip over, or properly analyze the pedestal events
- run->SetEventRange(1, MaxEvent); // Physics Event number, does not include scaler or control events.
- run->SetNscan(1);
- run->SetDataRequired(0x7);
- run->Print();
-
- // Define the analysis parameters
- TString ROOTFileName = Form(ROOTFileNamePattern, RunNumber, MaxEvent);
- analyzer->SetCountMode(2); // 0 = counter is # of physics triggers
- // 1 = counter is # of all decode reads
- // 2 = counter is event number
-
- analyzer->SetEvent(event);
- // Set EPICS event type
- analyzer->SetEpicsEvtType(180);
- // Define crate map
- analyzer->SetCrateMapFileName("MAPS/db_cratemap.dat");
- // Define output ROOT file
- analyzer->SetOutFile(ROOTFileName.Data());
- // Define DEF-file+
- analyzer->SetOdefFile("DEF-files/COIN/PRODUCTION/coin_production_hElec_pProt.def");
- // Define cuts file
- analyzer->SetCutFile("DEF-files/COIN/PRODUCTION/CUTS/coin_production_cuts.def"); // optional
- // File to record accounting information for cuts
- analyzer->SetSummaryFile(Form("REPORT_OUTPUT/COIN/PRODUCTION/summary_production_%d_%d.report", RunNumber, MaxEvent)); // optional
- // Start the actual analysis.
- analyzer->Process(run);
- // Create report file from template
- analyzer->PrintReport("TEMPLATES/COIN/PRODUCTION/coin_production.template",
- Form("REPORT_OUTPUT/COIN/PRODUCTION/replay_coin_production_%d_%d.report", RunNumber, MaxEvent)); // optional
-
-}
diff --git a/tests/elastic_test.cxx b/tests/elastic_test.cxx
deleted file mode 100644
index 4c9ad79f09c3cfefec9d399f8fd68d6718ee3551..0000000000000000000000000000000000000000
--- a/tests/elastic_test.cxx
+++ /dev/null
@@ -1,446 +0,0 @@
-#include "nlohmann/json.hpp"
-#include <cmath>
-#include <iostream>
-
-#include "ROOT/RDataFrame.hxx"
-#include "ROOT/RVec.hxx"
-
-#include "Math/Vector3D.h"
-#include "Math/Vector4D.h"
-#include "Math/VectorUtil.h"
-#include "TCanvas.h"
-#include "TLatex.h"
-#include "TStyle.h"
-#include "TSystem.h"
-R__LOAD_LIBRARY(libMathMore.so)
-R__LOAD_LIBRARY(libGenVector.so)
-
-R__LOAD_LIBRARY(libfmt.so)
-#include "fmt/core.h"
-
-#include "THStack.h"
-
-#ifdef __cpp_lib_filesystem
-#include <filesystem>
-namespace fs = std::filesystem;
-#else
-#include <experimental/filesystem>
-namespace fs = std::experimental::filesystem;
-#endif
-
-using RDFNode = decltype(ROOT::RDataFrame{0}.Filter(""));
-using Histo1DProxy =
- decltype(ROOT::RDataFrame{0}.Histo1D(ROOT::RDF::TH1DModel{"", "", 128u, 0., 0.}, ""));
-
-struct RDFInfo {
- RDFNode& df;
- const std::string title;
- RDFInfo(RDFNode& df, std::string_view title) : df{df}, title{title} {}
-};
-
-constexpr const double M_P = .938272;
-constexpr const double M_P2 = M_P * M_P;
-constexpr const double M_pion = 0.139;
-constexpr const double M_pion2 = M_pion * M_pion;
-constexpr const double M_e = .000511;
-
-// =================================================================================
-// Cuts
-// =================================================================================
-std::string goodTrackSHMS = "P.gtr.dp > -10 && P.gtr.dp < 22";
-std::string goodTrackHMS = "H.gtr.dp > -8 && H.gtr.dp < 8";
-
-std::string piCutSHMS = "P.cal.etottracknorm<1.0";
-//std::string piCutSHMS = "P.aero.npeSum > 1.0 && P.cal.eprtracknorm < 0.2 && P.cal.etottracknorm<1.0";
-
-std::string eCutHMS = "H.cal.etottracknorm > 0.50 && H.cal.etottracknorm < 2. && "
- "H.cer.npeSum > 1.";
-
-std::string epiCut = "P.aero.npeSum > 1.0 && P.cal.eprtracknorm < 0.2 && "
- "H.cer.npeSum > 1.0 && H.cal.etottracknorm > 0.6 && "
- "H.cal.etottracknorm < 2.0 && P.cal.etottracknorm<1.0";
-
-using Pvec3D = ROOT::Math::XYZVector;
-using Pvec4D = ROOT::Math::PxPyPzMVector;
-
-// =================================================================================
-// reconstruction
-// =================================================================================
-auto p_pion = [](double px, double py, double pz) {
- return Pvec4D{px * 0.996, py * 0.996, pz * 0.996, M_e};
-};
-auto p_electron = [](double px, double py, double pz) {
- return Pvec4D{px * 0.994, py * 0.994, pz * 0.994, M_e};
-};
-auto t = [](const double Egamma, Pvec4D& jpsi) {
- Pvec4D beam{0, 0, Egamma, 0};
- return (beam - jpsi).M2();
-};
-
-bool root_file_exists(std::string rootfile) {
- bool found_good_file = false;
- if (!gSystem->AccessPathName(rootfile.c_str())) {
- TFile file(rootfile.c_str());
- if (file.IsZombie()) {
- std::cout << rootfile << " is a zombie.\n";
- std::cout
- << " Did your replay finish? Check that the it is done before running this script.\n";
- return false;
- // return;
- } else {
- std::cout << " using : " << rootfile << "\n";
- return true;
- }
- }
- return false;
-}
-
-void elastic_test(int RunNumber = 6012, int nevents = 50000, int prompt = 0, int update = 0,
- int default_count_goal = 10000, int redo_timing = 0) {
-
- // ===============================================================================================
- // Initialization
- // ===============================================================================================
- using nlohmann::json;
- json j;
- {
- std::ifstream json_input_file("db2/run_list_coin.json");
- try {
- json_input_file >> j;
- } catch (json::parse_error) {
- std::cerr << "error: json file, db2/run_list.json, is incomplete or has broken syntax.\n";
- std::quick_exit(-127);
- }
- }
-
- auto runnum_str = std::to_string(RunNumber);
- if (j.find(runnum_str) == j.end()) {
- std::cout << "Run " << RunNumber << " not found in db2/run_list_coin.json\n";
- std::cout << "Check that run number and replay exists. \n";
- std::cout << "If problem persists please contact Sylvester (217-848-0565)\n";
- }
- double P0_shms_setting = j[runnum_str]["spectrometers"]["shms_momentum"].get<double>();
- double P0_shms = std::abs(P0_shms_setting);
-
- bool found_good_file = false;
-
- std::string rootfile =
- fmt::format("full_online/coin_replay_production_{}_{}.root", RunNumber, nevents);
- found_good_file = root_file_exists(rootfile.c_str());
- if (!found_good_file) {
- rootfile =
- fmt::format("ROOTfiles_volatile/coin_replay_production_{}_{}.root", RunNumber, nevents);
- found_good_file = root_file_exists(rootfile.c_str());
- }
- if (!found_good_file) {
- rootfile = fmt::format("ROOTfiles_csv/coin_replay_production_{}_{}.root", RunNumber, nevents);
- found_good_file = root_file_exists(rootfile.c_str());
- }
- if (!found_good_file) {
- rootfile = fmt::format("ROOTfiles/coin_replay_production_{}_{}.root", RunNumber, nevents);
- found_good_file = root_file_exists(rootfile.c_str());
- }
- if (!found_good_file) {
- std::cout << " Error: suitable root file not found\n";
- return;
- }
-
- // ===============================================================================================
- // Dataframe
- // ===============================================================================================
-
- ROOT::EnableImplicitMT(24);
-
- //---------------------------------------------------------------------------
- // Detector tree
- ROOT::RDataFrame d("T", rootfile);
-
- //// SHMS Scaler tree
- //ROOT::RDataFrame d_sh("TSP", rootfile);
- //// int N_scaler_events = *(d_sh.Count());
-
- auto d_coin = d.Filter("fEvtHdr.fEvtType == 4");
-
- // Good track cuts
- auto dHMSGoodTrack = d_coin.Filter(goodTrackHMS);
- auto dSHMSGoodTrack = d_coin.Filter(goodTrackSHMS);
- auto dCOINGoodTrack = dHMSGoodTrack.Filter(goodTrackSHMS)
- .Define("p_electron", p_electron, {"H.gtr.px", "H.gtr.py", "H.gtr.pz"})
- .Define("p_pion", p_pion, {"P.gtr.px", "P.gtr.py", "P.gtr.pz"});
- // PID cuts
- auto dHMSEl = dHMSGoodTrack.Filter(eCutHMS);
- auto dSHMSEl = dSHMSGoodTrack.Filter(piCutSHMS);
- auto dCOINEl = dCOINGoodTrack.Filter(eCutHMS + " && " + piCutSHMS);
- //.Filter(
- // [=](double npe, double dp) {
- // double p_track = P0_shms * (100.0 + dp) / 100.0;
- // // no cerenkov cut needed when momentum is below 2.8 GeV/c
- // if (p_track < 2.8) {
- // return true;
- // }
- // return npe > 1.0;
- // },
- // {"P.hgcer.npeSum", "P.gtr.dp"});
-
- // Timing cuts
- // Find the timing peak
- // Find the coin peak
- double coin_peak_center = 0;
- if (redo_timing) {
- auto h_coin_time =
- dCOINEl.Histo1D({"coin_time", "coin_time", 8000, 0, 1000}, "CTime.ePositronCoinTime_ROC2");
- h_coin_time->DrawClone();
- int coin_peak_bin = h_coin_time->GetMaximumBin();
- coin_peak_center = h_coin_time->GetBinCenter(coin_peak_bin);
- std::cout << "COINCIDENCE time peak found at: " << coin_peak_center << std::endl;
- } else {
- //coin_peak_center = 43.0; // run 7240-7241
- coin_peak_center = 23.0; // run 6012
- std::cout << "COINCIDENCE time peak: using pre-calculated value at: " << coin_peak_center
- << std::endl;
- ;
- }
- // timing cut lambdas
- // TODO: evaluate timing cut and offset for random background
- auto timing_cut = [=](double coin_time) { return std::abs(coin_time - coin_peak_center) < 2.; };
- // anti-timing set to 5x width of regular
- auto anti_timing_cut = [=](double coin_time) {
- return std::abs(coin_time - coin_peak_center - 28.) < 10.;
- };
-
- // timing counts
- auto dHMSElInTime = dHMSEl.Filter(timing_cut, {"CTime.ePositronCoinTime_ROC2"});
- auto dHMSElRandom = dHMSEl.Filter(anti_timing_cut, {"CTime.ePositronCoinTime_ROC2"});
- auto dSHMSElInTime = dSHMSEl.Filter(timing_cut, {"CTime.ePositronCoinTime_ROC2"});
- auto dSHMSElRandom = dSHMSEl.Filter(anti_timing_cut, {"CTime.ePositronCoinTime_ROC2"});
-
- auto dCOINElInTime = dCOINEl.Filter(timing_cut, {"CTime.ePiCoinTime_ROC2"});
- auto dCOINElRandom = dCOINEl.Filter(anti_timing_cut, {"CTime.ePiCoinTime_ROC2"});
-
- // Output root file
- //auto out_file =
- // new TFile(fmt::format("monitoring/{}/good_csv_counter.root", RunNumber).c_str(), "UPDATE");
- //out_file->cd();
-
- // =========================================================================================
- // Histograms
- // =========================================================================================
- // 2D correlations
- auto hTheta2DNoCuts = d_coin.Histo2D(
- {"theta2D", "No cuts;#theta_{SHMS};#theta_{HMS};#counts", 50, -.1, .1, 50, -.1, .1},
- "P.gtr.th", "H.gtr.th");
- auto hTheta2DTracking = dCOINGoodTrack.Histo2D(
- {"theta2D", "Cuts: tracking;#theta_{SHMS};#theta_{HMS};#counts", 50, -.1, .1, 50, -.1, .1},
- "P.gtr.th", "H.gtr.th");
- auto hTheta2DPID =
- dCOINEl.Histo2D({"theta2D", "Cuts: tracking+PID;#theta_{SHMS};#theta_{HMS};#counts", 50, -.1,
- .1, 50, -.1, .1},
- "P.gtr.th", "H.gtr.th");
- auto hTheta2DTiming =
- dCOINElInTime.Histo2D({"theta2D", "Cuts: tracking+PID;#theta_{SHMS};#theta_{HMS};#counts", 50,
- -.1, .1, 50, -.1, .1},
- "P.gtr.th", "H.gtr.th");
- // timing
- auto hCoinTimeNoCuts =
- d_coin.Histo1D({"coin_time.NoCuts", "No Cuts;coin_time;counts", 8000, 0, 1000},
- "CTime.ePositronCoinTime_ROC2");
- auto hCoinTimeTracking = dCOINGoodTrack.Histo1D(
- {"coin_time.Tracking", "Cuts: Tracking;coin_time;counts", 8000, 0, 1000},
- "CTime.ePositronCoinTime_ROC2");
- auto hCoinTimePID =
- dCOINEl.Histo1D({"coin_time.PID", "Cuts: Tracking+PID;coin_time;counts", 8000, 0, 1000},
- "CTime.ePositronCoinTime_ROC2");
- auto hCoinTimeTiming = dCOINElInTime.Histo1D(
- {"coin_time.Timing", "Cuts: Tracking+PID+Timing;coin_time;counts", 8000, 0, 1000},
- "CTime.ePositronCoinTime_ROC2");
-
- auto hRandCoinTimePID = dCOINElRandom.Histo1D(
- {"rand_coin_time.PID", "Cuts: Tracking+PID;coin_time;counts", 8000, 0, 1000},
- "CTime.ePositronCoinTime_ROC2");
-
- // P.gtr.dp
- auto hPdpNoCuts =
- d_coin.Histo1D({"P.gtr.dp.NoCuts", "No Cuts;#deltap [%];counts", 200, -30, 40}, "P.gtr.dp");
- auto hPdpTracking = dSHMSGoodTrack.Histo1D(
- {"P.gtr.dp.Tracking", "Cuts: Tracking;#deltap [%];counts", 200, -30, 40}, "P.gtr.dp");
- auto hPdpPID = dSHMSEl.Histo1D(
- {"P.gtr.dp.PID", "Cuts: Tracking+PID;#deltap [%];counts", 200, -30, 40}, "P.gtr.dp");
- auto hPdpTiming = dSHMSElInTime.Histo1D(
- {"P.gtr.dp.Timing", "Cuts: Tracking+PID+Timing;#deltap [%];counts", 200, -30, 40},
- "P.gtr.dp");
- // P.gtr.th
- auto hPthNoCuts = d_coin.Histo1D(
- {"P.gtr.th.NoCuts", "No Cuts;#theta_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.th");
- auto hPthTracking = dSHMSGoodTrack.Histo1D(
- {"P.gtr.th.Tracking", "Cuts: Tracking;#theta_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.th");
- auto hPthPID = dSHMSEl.Histo1D(
- {"P.gtr.th.PID", "Cuts: Tracking+PID;#theta_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.th");
- auto hPthTiming = dSHMSElInTime.Histo1D(
- {"P.gtr.th.Timing", "Cuts: Tracking+PID+Timing;#theta_{SHMS};counts", 200, -0.1, 0.1},
- "P.gtr.th");
- // P.gtr.ph
- auto hPphNoCuts =
- d_coin.Histo1D({"P.gtr.ph.NoCuts", "No Cuts;#phi_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.ph");
- auto hPphTracking = dSHMSGoodTrack.Histo1D(
- {"P.gtr.ph.Tracking", "Cuts: Tracking;#phi_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.ph");
- auto hPphPID = dSHMSEl.Histo1D(
- {"P.gtr.ph.PID", "Cuts: Tracking+PID;#phi_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.ph");
- auto hPphTiming = dSHMSElInTime.Histo1D(
- {"P.gtr.ph.Timing", "Cuts: Tracking+PID+Timing;#phi_{SHMS};counts", 200, -0.1, 0.1},
- "P.gtr.ph");
- // P.gtr.y
- auto hPyNoCuts =
- d_coin.Histo1D({"P.gtr.y.NoCuts", "No Cuts;ytar;counts", 200, -10., 10.}, "P.gtr.y");
- auto hPyTracking = dSHMSGoodTrack.Histo1D(
- {"P.gtr.y.Tracking", "Cuts: Tracking;ytar;counts", 200, -10., 10.}, "P.gtr.y");
- auto hPyPID =
- dSHMSEl.Histo1D({"P.gtr.y.PID", "Cuts: Tracking+PID;ytar;counts", 200, -10., 10.}, "P.gtr.y");
- auto hPyTiming = dSHMSElInTime.Histo1D(
- {"P.gtr.y.Timing", "Cuts: Tracking+PID+Timing;ytar;counts", 200, -10., 10.}, "P.gtr.y");
- // P.cal.etottracknorm
- auto hPcalEPNoCuts =
- d_coin.Histo1D({"P.cal.etottracknorm.NoCuts", "No Cuts;SHMS E/P;counts", 200, -.5, 1.5},
- "P.cal.etottracknorm");
- auto hPcalEPTracking = dSHMSGoodTrack.Histo1D(
- {"P.cal.etottracknorm.Tracking", "Cuts: Tracking;SHMS E/P;counts", 200, -.5, 1.5},
- "P.cal.etottracknorm");
- auto hPcalEPPID = dSHMSEl.Histo1D(
- {"P.cal.etottracknorm.PID", "Cuts: Tracking+PID;SHMS E/P;counts", 200, -.5, 1.5},
- "P.cal.etottracknorm");
- auto hPcalEPAll = dCOINElInTime.Histo1D(
- {"P.cal.etottracknorm.All", "Cuts: Tracking+PID+Coincidence;SHMS E/P;counts", 200, -.5, 1.5},
- "P.cal.etottracknorm");
- // P.ngcer.npeSum
- auto hPcerNpheNoCuts = d_coin.Histo1D(
- {"P.ngcer.npeSum.NoCuts", "No Cuts;SHMS NGC #phe;counts", 200, -5, 76}, "P.ngcer.npeSum");
- auto hPcerNpheTracking = dSHMSGoodTrack.Histo1D(
- {"P.ngcer.npeSum.Tracking", "Cuts: Tracking;SHMS NGC #phe;counts", 200, -5, 76},
- "P.ngcer.npeSum");
- auto hPcerNphePID = dSHMSEl.Histo1D(
- {"P.ngcer.npeSum.PID", "Cuts: Tracking+PID;SHMS NGC #phe;counts", 200, -5, 76},
- "P.ngcer.npeSum");
- auto hPcerNpheAll = dCOINElInTime.Histo1D(
- {"P.ngcer.npeSum.All", "Cuts: Tracking+PID+Coincidence;SHMS NGC #phe;counts", 200, -5, 76},
- "P.ngcer.npeSum");
- // P.hgcer.npeSum
- auto hPhgcerNpheNoCuts = d_coin.Histo1D(
- {"P.hgcer.npeSum.NoCuts", "No Cuts;SHMS HGC #phe;counts", 200, -5, 76}, "P.hgcer.npeSum");
- auto hPhgcerNpheTracking = dSHMSGoodTrack.Histo1D(
- {"P.hgcer.npeSum.Tracking", "Cuts: Tracking;SHMS HGC #phe;counts", 200, -5, 76},
- "P.hgcer.npeSum");
- auto hPhgcerNphePID = dSHMSEl.Histo1D(
- {"P.hgcer.npeSum.PID", "Cuts: Tracking+PID;SHMS HGC #phe;counts", 200, -5, 76},
- "P.hgcer.npeSum");
- auto hPhgcerNpheAll = dCOINElInTime.Histo1D(
- {"P.hgcer.npeSum.All", "Cuts: Tracking+PID+Coincidence;SHMS HGC #phe;counts", 200, -5, 76},
- "P.hgcer.npeSum");
- // H.cal.etottracknorm
- auto hHcalEPNoCuts =
- d_coin.Histo1D({"H.cal.etottracknorm.NoCuts", "No Cuts;HMS E/P;counts", 200, -.5, 1.5},
- "H.cal.etottracknorm");
- auto hHcalEPTracking = dHMSGoodTrack.Histo1D(
- {"H.cal.etottracknorm.Tracking", "Cuts: Tracking;HMS E/P;counts", 200, -.5, 1.5},
- "H.cal.etottracknorm");
- auto hHcalEPPID = dHMSEl.Histo1D(
- {"H.cal.etottracknorm.PID", "Cuts: Tracking+PID;HMS E/P;counts", 200, -.5, 1.5},
- "H.cal.etottracknorm");
- auto hHcalEPAll = dCOINElInTime.Histo1D(
- {"H.cal.etottracknorm.All", "Cuts: Tracking+PID+Coincidence;HMS E/P;counts", 200, -.5, 1.5},
- "H.cal.etottracknorm");
- // H.cer.npeSum
- auto hHcerNpheNoCuts = d_coin.Histo1D(
- {"H.cer.npeSum.NoCuts", "No Cuts;HMS Cer #phe;counts", 200, -1, 15}, "H.cer.npeSum");
- auto hHcerNpheTracking = dHMSGoodTrack.Histo1D(
- {"H.cer.npeSum.Tracking", "Cuts: Tracking;HMS Cer #phe;counts", 200, -1, 15}, "H.cer.npeSum");
- auto hHcerNphePID = dHMSEl.Histo1D(
- {"H.cer.npeSum.PID", "Cuts: Tracking+PID+Coincidence;HMS Cer #phe;counts", 200, -1, 15},
- "H.cer.npeSum");
- auto hHcerNpheAll = dCOINElInTime.Histo1D(
- {"H.cer.npeSum.PID", "Cuts: Tracking+PID+Coincidence;HMS Cer #phe;counts", 200, -1, 15},
- "H.cer.npeSum");
- // H.gtr.dp
- auto hHdpNoCuts =
- d_coin.Histo1D({"H.gtr.dp.NoCuts", "No Cuts;#deltap [%];counts", 200, -30, 40}, "H.gtr.dp");
- auto hHdpTracking = dHMSGoodTrack.Histo1D(
- {"H.gtr.dp.Tracking", "Cuts: Tracking;#deltap [%];counts", 200, -30, 40}, "H.gtr.dp");
- auto hHdpPID = dHMSEl.Histo1D(
- {"H.gtr.dp.PID", "Cuts: Tracking+PID;#deltap [%];counts", 200, -30, 40}, "H.gtr.dp");
- auto hHdpTiming = dHMSElInTime.Histo1D(
- {"H.gtr.dp.Timing", "Cuts: Tracking+PID+Timing;#deltap [%];counts", 200, -30, 40},
- "H.gtr.dp");
- // H.gtr.th
- auto hHthNoCuts = d_coin.Histo1D(
- {"H.gtr.th.NoCuts", "No Cuts;#theta_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.th");
- auto hHthTracking = dHMSGoodTrack.Histo1D(
- {"H.gtr.th.Tracking", "Cuts: Tracking;#theta_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.th");
- auto hHthPID = dHMSEl.Histo1D(
- {"H.gtr.th.PID", "Cuts: Tracking+PID;#theta_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.th");
- auto hHthTiming = dHMSElInTime.Histo1D(
- {"H.gtr.th.Timing", "Cuts: Tracking+PID+Timing;#theta_{HMS};counts", 200, -0.1, 0.1},
- "H.gtr.th");
- // H.gtr.ph
- auto hHphNoCuts =
- d_coin.Histo1D({"H.gtr.ph.NoCuts", "No Cuts;#phi_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.ph");
- auto hHphTracking = dHMSGoodTrack.Histo1D(
- {"H.gtr.ph.Tracking", "Cuts: Tracking;#phi_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.ph");
- auto hHphPID = dHMSEl.Histo1D(
- {"H.gtr.ph.PID", "Cuts: Tracking+PID;#phi_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.ph");
- auto hHphTiming = dHMSElInTime.Histo1D(
- {"H.gtr.ph.Timing", "Cuts: Tracking+PID+Timing;#phi_{HMS};counts", 200, -0.1, 0.1},
- "H.gtr.ph");
- // H.gtr.y
- auto hHyNoCuts =
- d_coin.Histo1D({"H.gtr.y.NoCuts", "No Cuts;ytar;counts", 200, -10., 10.}, "H.gtr.y");
- auto hHyTracking = dHMSGoodTrack.Histo1D(
- {"H.gtr.y.Tracking", "Cuts: Tracking;ytar;counts", 200, -10., 10.}, "H.gtr.y");
- auto hHyPID =
- dHMSEl.Histo1D({"H.gtr.y.PID", "Cuts: Tracking+PID;ytar;counts", 200, -10., 10.}, "H.gtr.y");
- auto hHyTiming = dHMSElInTime.Histo1D(
- {"H.gtr.y.Timing", "Cuts: Tracking+PID+Timing;ytar;counts", 200, -10., 10.}, "H.gtr.y");
-
- // scalers
- //auto total_charge = d_sh.Max("P.BCM4B.scalerChargeCut");
- //auto shms_el_real_scaler = d_sh.Max("P.pEL_REAL.scaler");
- //auto hms_el_real_scaler = d_sh.Max("P.hEL_REAL.scaler");
- //auto time_1MHz = d_sh.Max("P.1MHz.scalerTime");
- //auto time_1MHz_cut = d_sh.Max("P.1MHz.scalerTimeCut");
-
- auto yield_all = d.Count();
- // 5 timing cut widths worth of random backgrounds
- auto yield_coin = d_coin.Count();
- auto yield_HMSGoodTrack = dHMSGoodTrack.Count();
- auto yield_SHMSGoodTrack = dSHMSGoodTrack.Count();
- auto yield_COINGoodTrack = dCOINGoodTrack.Count();
- auto yield_HMSEl = dHMSEl.Count();
- auto yield_SHMSEl = dSHMSEl.Count();
- auto yield_COINEl = dCOINEl.Count();
- auto yield_HMSElInTime = dHMSElInTime.Count();
- auto yield_HMSElRandom = dHMSElRandom.Count();
- auto yield_SHMSElInTime = dSHMSElInTime.Count();
- auto yield_SHMSElRandom = dSHMSElRandom.Count();
- auto yield_COINElInTime = dCOINElInTime.Count();
- auto yield_COINElRandom = dCOINElRandom.Count();
- auto yield_coin_raw = dCOINElInTime.Count();
- auto yield_coin_random = dCOINElRandom.Count();
-
-
- // -------------------------------------
- // End lazy eval
- // -------------------------------------
- auto n_coin_good = *yield_coin_raw - *yield_coin_random / 5.;
- auto n_HMSElGood = *yield_HMSElInTime - *yield_HMSElRandom / 5;
- auto n_SHMSElGood = *yield_SHMSElInTime - *yield_SHMSElRandom / 5;
- auto n_COINElGood = *yield_COINElInTime - *yield_COINElRandom / 5;
-
- hPdpNoCuts->DrawCopy();
- //std::cout << " coin COUNTS : " << *(d_coin.Count()) << "\n";
- //std::cout << " yield_HMSEl : " << *(yield_HMSEl) << "\n";
- std::cout << " yield_COINEl : " << *(yield_COINEl) << "\n";
- //std::cout << " ALL COUNTS : " << *yield_all << "\n";
- //std::cout << " GOOD COUNTS : " << n_COINElGood << "\n";
- //
- if( 4 != (*(yield_COINEl)) ){
- std::exit(-1);
- }
- std::exit(0);
-}
diff --git a/tests/elastic_test.sh b/tests/elastic_test.sh
deleted file mode 100644
index 5d137e59c5413ab3ab02a67612ca8466182b04f9..0000000000000000000000000000000000000000
--- a/tests/elastic_test.sh
+++ /dev/null
@@ -1,26 +0,0 @@
-#!/bin/bash
-
-#echo "This is the elastic testing..."
-#echo " "
-#echo "There are currently 0 tests to run!"
-#which hcana
-#
-#ls -lrth
-#ls -lrth build
-#
-#git clone git@eicweb.phy.anl.gov:jlab/hallc/exp/CSV/hallc_replay_csv.git
-#git clone git@eicweb.phy.anl.gov:jlab/hallc/exp/CSV/online_csv.git
-#cd online_csv
-#ln -s ../hallc_reaply_csv/PARAM
-## and the reset
-#
-#mkdir raw
-#pushd raw
-# wget coin.dat
-#popd
-
-
-singularity help build/Singularity.hcana.simg
-
-
-singularity exec build/Singularity.hcana.simg hcana tests/elastic_test.cxx
diff --git a/tests/elastic_test2.sh b/tests/elastic_test2.sh
deleted file mode 100644
index 36da2058d88749b20c9af12ce265f23befd157c2..0000000000000000000000000000000000000000
--- a/tests/elastic_test2.sh
+++ /dev/null
@@ -1,30 +0,0 @@
-#!/bin/bash
-
-echo "This is the elastic testing..."
-echo " "
-echo "There are currently 0 tests to run!"
-which hcana
-
-ls -lrth
-ls -lrth build
-
-git clone git@eicweb.phy.anl.gov:jlab/hallc/exp/CSV/hallc_replay_csv.git
-git clone git@eicweb.phy.anl.gov:jlab/hallc/exp/CSV/online_csv.git
-
-cd online_csv
-ln -s ../hallc_reaply_csv/PARAM
-# and the reset
-
-mkdir raw
-pushd raw
- wget coin.dat
-popd
-
-
-singularity help build/Singularity.hcana.simg
-
-singularity exec build/Singularity.hcana.simg which hcana
-
-singularity exec build/Singularity.hcana.simg hcana tests/my_root_script.cxx
-
-echo " WOOOO"
diff --git a/tests/my_root_script.cxx b/tests/my_root_script.cxx
deleted file mode 100644
index 0d3686284d799790b779b7724cdbca81d3a5c129..0000000000000000000000000000000000000000
--- a/tests/my_root_script.cxx
+++ /dev/null
@@ -1,20 +0,0 @@
-void my_root_script() {
-
-
-
- std::cout << "Hello from my_root_script.cxx!\n";
-
- std::cout << "This should be run with singularity\n";
- double pi = 3.14;
-
- auto pi_equals_3 = (3 == pi);
- std::cout << " pi_equals_3 = " << pi_equals_3 << "\n";
-
- if( pi_equals_3) {
- std::cout << "what the hell?\n";
- std::exit( 0 );
- }
- /* else */
-
- std::exit( -1 );
-}
diff --git a/tests/replay_elastic_data.sh b/tests/replay_elastic_data.sh
deleted file mode 100644
index 48b147f737be5ce6b44ed32e39f3505b377dccb8..0000000000000000000000000000000000000000
--- a/tests/replay_elastic_data.sh
+++ /dev/null
@@ -1,38 +0,0 @@
-#!/bin/bash
-
-echo "This is the elastic testing..."
-echo " "
-echo "There are currently 0 tests to run!"
-which hcana
-
-ls -lrth
-ls -lrth build
-mkdir -p ROOTfiles
-git clone https://eicweb.phy.anl.gov/whit/ci_test_data.git
-git clone https://eicweb.phy.anl.gov/jlab/hallc/exp/CSV/hallc_replay_csv.git
-git clone https://eicweb.phy.anl.gov/jlab/hallc/exp/CSV/online_csv.git
-
-cd online_csv
-mkdir -p logs
-ln -s ../hallc_replay_csv/PARAM
-ln -s ../hallc_replay_csv/DBASE
-ln -s ../hallc_replay_csv/CALIBRATION
-ln -s ../hallc_replay_csv/DEF-files
-ln -s ../hallc_replay_csv/MAPS
-ln -s ../hallc_replay_csv/SCRIPTS
-ln -s ../hallc_replay_csv/DATFILES
-ln -s ../ci_test_data/raw
-ln -s ../ROOTfiles
-# and the reset
-
-ls -lrth
-ls -lrth raw/
-ls -lrth ROOTfiles/
-pwd
-# run replay script
-df -h
-
-singularity exec ../build/Singularity.hcana.simg hcana -b -q "../tests/elastic_coin_replay.cxx+(6012,50000)"
-
-echo "hcana calls... the coin replay script and outputs blah.root"
-