diff --git a/dvcs/config.yml b/dvcs/config.yml
index 0bd08d37a8510efd65e100c5b57d76a311bf415a..41b2ee757d3765c85def0ef0bbd20f53821f6624 100644
--- a/dvcs/config.yml
+++ b/dvcs/config.yml
@@ -12,10 +12,11 @@ dvcs:results:
stage: collect
needs: ["dvcs:process"]
script:
- - echo "All DVCS benchmarks successful"
+ - pip install junitparser
+ - python dvcs/scripts/merge_results.py
artifacts:
paths:
- results
reports:
- junit: ["results/dvcs/report.xml"]
+ junit: ["results/dvcs/dvcs_report.xml"]
diff --git a/dvcs/dvcs.sh b/dvcs/dvcs.sh
index ebab4c723d52217814fd939eec700c88319fa62a..436406f26997f0f265f342f5222b00c1486493f1 100644
--- a/dvcs/dvcs.sh
+++ b/dvcs/dvcs.sh
@@ -26,6 +26,8 @@ if [[ "$?" -ne "0" ]] ; then
exit 1
fi
+export JUGGLER_N_EVENTS=10
+
## run geant4 simulations
npsim --runType batch \
--part.minimalKineticEnergy 1000*GeV \
@@ -48,12 +50,12 @@ if [[ "$?" -ne "0" ]] ; then
fi
mkdir -p results/dvcs
-echo "STAND-IN FOR ANALYSIS SCRIPT"
-#root -b -q "dis/scripts/rec_dis_electrons.cxx(\"${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}\")"
-#if [[ "$?" -ne "0" ]] ; then
-# echo "ERROR running root script"
-# exit 1
-#fi
+
+root -b -q "dvcs/scripts/dvcs_tests.cxx(\"${JUGGLER_REC_FILE}\")"
+if [[ "$?" -ne "0" ]] ; then
+ echo "ERROR running root script"
+ exit 1
+fi
# copy data if it is not too big
if [[ "${JUGGLER_N_EVENTS}" -lt "500" ]] ; then
@@ -61,8 +63,8 @@ cp ${JUGGLER_REC_FILE} results/dvcs/.
fi
# Collect the results
-cp dvcs/report.xml results/dvcs/.
-cp dvcs/report2.xml results/dvcs/.
+#cp dvcs/report.xml results/dvcs/.
+#cp dvcs/report2.xml results/dvcs/.
diff --git a/dvcs/scripts/dvcs_tests.cxx b/dvcs/scripts/dvcs_tests.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..4f37494c040927f282f28193f45b5862c88f601c
--- /dev/null
+++ b/dvcs/scripts/dvcs_tests.cxx
@@ -0,0 +1,126 @@
+#include <ROOT/RDataFrame.hxx>
+#include <cmath>
+#include "fmt/color.h"
+#include "fmt/core.h"
+#include <iostream>
+#include <string>
+#include <vector>
+#include <nlohmann/json.hpp>
+using json = nlohmann::json;
+
+#include "dd4pod/Geant4ParticleCollection.h"
+#include "eicd/TrackParametersCollection.h"
+#include "eicd/ClusterCollection.h"
+#include "eicd/ClusterData.h"
+
+using ROOT::RDataFrame;
+using namespace ROOT::VecOps;
+
+auto p_track = [](std::vector<eic::TrackParametersData> const& in) {
+ std::vector<double> result;
+ for (size_t i = 0; i < in.size(); ++i) {
+ result.push_back(std::abs(1.0/(in[i].qOverP)));
+ }
+ return result;
+};
+
+
+auto pt = [](std::vector<dd4pod::Geant4ParticleData> const& in){
+ std::vector<float> result;
+ for (size_t i = 0; i < in.size(); ++i) {
+ result.push_back(std::sqrt(in[i].psx * in[i].psx + in[i].psy * in[i].psy));
+ }
+ return result;
+};
+
+auto momentum = [](std::vector<ROOT::Math::PxPyPzMVector> const& in) {
+ std::vector<double> result;
+ for (size_t i = 0; i < in.size(); ++i) {
+ result.push_back(in[i].E());
+ }
+ return result;
+};
+auto theta = [](std::vector<ROOT::Math::PxPyPzMVector> const& in) {
+ std::vector<double> result;
+ for (size_t i = 0; i < in.size(); ++i) {
+ result.push_back(in[i].Theta()*180/M_PI);
+ }
+ return result;
+};
+auto fourvec = [](ROOT::VecOps::RVec<dd4pod::Geant4ParticleData> const& in) {
+ std::vector<ROOT::Math::PxPyPzMVector> result;
+ ROOT::Math::PxPyPzMVector lv;
+ for (size_t i = 0; i < in.size(); ++i) {
+ lv.SetCoordinates(in[i].psx, in[i].psy, in[i].psz, in[i].mass);
+ result.push_back(lv);
+ }
+ return result;
+};
+
+auto delta_p = [](const std::vector<double>& tracks, const std::vector<double>& thrown) {
+ std::vector<double> res;
+ for (const auto& p1 : thrown) {
+ for (const auto& p2 : tracks) {
+ res.push_back(p1 - p2);
+ }
+ }
+ return res;
+};
+
+
+void dvcs_tests(const char* fname = "rec_dvcs.root"){
+
+ fmt::print(fmt::emphasis::bold | fg(fmt::color::forest_green), "Running DVCS analysis...\n");
+
+ // Run this in multi-threaded mode if desired
+ ROOT::EnableImplicitMT();
+ ROOT::RDataFrame df("events", fname);
+
+ using ROOT::Math::PxPyPzMVector;
+ PxPyPzMVector p_ebeam = {0,0,-10, 0.000511};
+ PxPyPzMVector p_pbeam = {0,0,275, 0.938 };
+
+ auto eprime = [](ROOT::VecOps::RVec<dd4pod::Geant4ParticleData> const& in) {
+ for(const auto& p : in){
+ if(p.pdgID == 11 ) {
+ return PxPyPzMVector(p.psx,p.psy,p.psz,p.mass);
+ }
+ }
+ return PxPyPzMVector(0,0,0,0);
+ };
+ auto q_vec = [=](PxPyPzMVector const& p) {
+ return p_ebeam - p;
+ };
+
+ auto df0 = df.Define("isThrown", "mcparticles2.genStatus == 1")
+ .Define("thrownParticles", "mcparticles2[isThrown]")
+ .Define("thrownP", fourvec, {"thrownParticles"})
+ .Define("p_thrown", momentum, {"thrownP"})
+ .Define("nTracks", "outputTrackParameters.size()")
+ .Define("p_track", p_track, {"outputTrackParameters"})
+ .Define("p_track1", p_track, {"outputTrackParameters1"})
+ .Define("p_track2", p_track, {"outputTrackParameters2"})
+ .Define("delta_p",delta_p, {"p_track", "p_thrown"})
+ .Define("delta_p1",delta_p, {"p_track1", "p_thrown"})
+ .Define("delta_p2",delta_p, {"p_track2", "p_thrown"})
+ .Define("eprime", eprime, {"thrownParticles"})
+ .Define("q", q_vec, {"eprime"})
+ .Define("Q2", "-1.0*(q.Dot(q))");
+
+ auto h_Q2 = df0.Histo1D({"h_Q2", "; Q^{2} [GeV^{2}/c^{2}]", 100, 0, 30}, "Q2");
+ auto n_Q2 = df0.Filter("Q2>1").Count();
+
+ auto c = new TCanvas();
+ h_Q2->DrawCopy();
+ c->SaveAs("results/dvcs/Q2.png");
+ c->SaveAs("results/dvcs/Q2.pdf");
+ fmt::print("{} DVCS events\n",*n_Q2);
+
+ // write output results to json file
+ json j;
+ j["Q2 cut"]["pass"] = *n_Q2;
+ j["Q2 cut"]["fail"] = 1;
+ std::ofstream o("results/dvcs/dvcs_tests.json");
+ o << std::setw(4) << j << std::endl;
+ //df0.Snapshot("testing","derp.root");
+}
diff --git a/dvcs/scripts/merge_results.py b/dvcs/scripts/merge_results.py
new file mode 100755
index 0000000000000000000000000000000000000000..3ec0cca730710239bcb9d674934a0d0dc0850d50
--- /dev/null
+++ b/dvcs/scripts/merge_results.py
@@ -0,0 +1,64 @@
+#!/usr/bin/env python3
+
+import json
+from junitparser import TestCase, TestSuite, JUnitXml, Skipped, Error, IntAttr, FloatAttr
+
+
+# Create the new element by subclassing Element or one of its child class,
+## and add custom attributes to it.
+#class MyTestCase(TestCase):
+# foo = Attr()
+
+# Add the custom attribute
+#TestCase.id = IntAttr('id')
+TestCase.efficiency = FloatAttr('efficiency')
+#TestCase.custom = Attr('custom')
+#case = TestCase()
+#case.id = 123
+#case.rate = 0.95
+#case.custom = 'foobar'
+
+# After looking at two different python libraries (junit-xml and junitparser)
+# junitparser looks the most robust
+# https://github.com/weiwei/junitparser
+
+def merge_results():
+ results = None;
+ with open("results/dvcs/dvcs_tests.json","r") as f:
+ results = json.load(f)
+
+ # Create suite and add cases
+ suite = TestSuite('dvcs')
+ suite.add_property('energy', '10-on-100')
+
+ for tname,tres in results.items():
+ for ttype, tval in tres.items():
+ # Create cases
+ case1 = TestCase(tname)
+ case1.time = 1.0
+ case1.efficiency = tval
+ case1.classname = ttype
+ suite.add_testcase(case1)
+
+ xml = JUnitXml()
+ xml.add_testsuite(suite)
+ xml.write('results/dvcs/dvcs_report.xml',pretty=True)
+
+
+#test code for junit-xml:
+#from junit_xml import TestSuite, TestCase
+# test_cases = []
+# print(test_name)
+# print(test_res)
+# for test_type, test_val in test_res.items():
+# test_cases.append(TestCase(test_name, "dvcs.dvcs_tests.{}".format(test_type), 10, str(test_val), 'I am stderr!'))
+# ts = TestSuite("my test suite", test_cases)
+# # pretty printing is on by default but can be disabled using prettyprint=False
+# print(TestSuite.to_xml_string([ts]))
+# # you can also write the XML to a file and not pretty print it
+# with open('results/dvcs/dvcs_report.xml', 'w') as f:
+# TestSuite.to_file(f, [ts], prettyprint=True)
+
+if __name__ == "__main__":
+ # execute only if run as a script
+ merge_results()