From 411559acdc2f154e6e2a0a6d279f21eb56a64a61 Mon Sep 17 00:00:00 2001
From: Dave Casper <dcasper@uci.edu>
Date: Wed, 8 Dec 2021 13:23:51 -0800
Subject: [PATCH] Fix integer overflow in persistent classes
---
.../src/CaloHitCollectionCnv.cxx | 8 +-
.../src/CaloHitCollectionCnv.h | 6 +-
.../CaloHits/CaloHitCollectionCnv_p1a.h | 41 +++
.../CaloHits/CaloHitCollection_p1a.h | 59 +++
.../CaloSimEventTPCnvDict.h | 2 +
.../FaserCaloSimEventTPCnv/selection.xml | 2 +
.../src/CaloHits/CaloHitCollectionCnv_p1a.cxx | 339 ++++++++++++++++++
.../src/NeutrinoHitCollectionCnv.cxx | 9 +-
.../src/NeutrinoHitCollectionCnv.h | 8 +-
.../NeutrinoHitCollectionCnv_p1a.h | 41 +++
.../NeutrinoHits/NeutrinoHitCollection_p1a.h | 57 +++
.../NeutrinoSimEventTPCnvDict.h | 3 +
.../NeutrinoSimEventTPCnv/selection.xml | 1 +
.../NeutrinoHitCollectionCnv_p1a.cxx | 334 +++++++++++++++++
14 files changed, 900 insertions(+), 10 deletions(-)
create mode 100644 Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1a.h
create mode 100644 Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1a.h
create mode 100644 Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/src/CaloHits/CaloHitCollectionCnv_p1a.cxx
create mode 100644 Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h
create mode 100644 Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1a.h
create mode 100644 Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/src/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.cxx
diff --git a/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.cxx b/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.cxx
index 2ac860d86..6abec9240 100644
--- a/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.cxx
+++ b/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.cxx
@@ -18,12 +18,16 @@ CaloHitCollection_PERS* CaloHitCollectionCnv::createPersistent(CaloHitCollection
CaloHitCollection* CaloHitCollectionCnv::createTransient() {
MsgStream mlog(msgSvc(), "CaloHitCollectionConverter" );
CaloHitCollectionCnv_p1 converter_p1;
+ CaloHitCollectionCnv_p1a converter_p1a;
static const pool::Guid p1_guid("134E8045-AB99-43EF-9AD1-324C48830B64");
- // static const pool::Guid p1_guid("B2573A16-4B46-4E1E-98E3-F93421680779");
+ static const pool::Guid p1a_guid("02C108EE-0AD9-4444-83BD-D5273FCDEF6F");
CaloHitCollection *trans_cont(0);
- if( this->compareClassGuid(p1_guid)) {
+ if( this->compareClassGuid(p1a_guid)) {
+ std::unique_ptr< CaloHitCollection_p1a > col_vect( this->poolReadObject< CaloHitCollection_p1a >() );
+ trans_cont = converter_p1a.createTransient( col_vect.get(), mlog );
+ } else if( this->compareClassGuid(p1_guid)) {
std::unique_ptr< CaloHitCollection_p1 > col_vect( this->poolReadObject< CaloHitCollection_p1 >() );
trans_cont = converter_p1.createTransient( col_vect.get(), mlog );
} else {
diff --git a/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.h b/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.h
index b46894c5f..1de7cc3e0 100644
--- a/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.h
+++ b/Calorimeter/CaloEventCnv/FaserCaloSimEventAthenaPool/src/CaloHitCollectionCnv.h
@@ -8,12 +8,14 @@
#include "FaserCaloSimEvent/CaloHitCollection.h"
#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1.h"
#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1.h"
+#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1a.h"
+#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1a.h"
#include "AthenaPoolCnvSvc/T_AthenaPoolCustomCnv.h"
// Gaudi
#include "GaudiKernel/MsgStream.h"
// typedef to the latest persistent version
-typedef CaloHitCollection_p1 CaloHitCollection_PERS;
-typedef CaloHitCollectionCnv_p1 CaloHitCollectionCnv_PERS;
+typedef CaloHitCollection_p1a CaloHitCollection_PERS;
+typedef CaloHitCollectionCnv_p1a CaloHitCollectionCnv_PERS;
class CaloHitCollectionCnv : public T_AthenaPoolCustomCnv<CaloHitCollection, CaloHitCollection_PERS > {
friend class CnvFactory<CaloHitCollectionCnv>;
diff --git a/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1a.h b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1a.h
new file mode 100644
index 000000000..80534aae7
--- /dev/null
+++ b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1a.h
@@ -0,0 +1,41 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef CALOHITCOLLECTIONCNV_P1A_H
+#define CALOHITCOLLECTIONCNV_P1A_H
+
+// CaloHitCollectionCnv_p1, T/P separation of Calo Hits
+// author D.Costanzo <davide.costanzo@cern.ch>
+// author O.Arnaez <olivier.arnaez@cern.ch>
+
+#include "AthenaPoolCnvSvc/T_AthenaPoolTPConverter.h"
+#include "FaserCaloSimEvent/CaloHitCollection.h"
+#include "CaloHitCollection_p1a.h"
+
+
+class CaloHitCollectionCnv_p1a : public T_AthenaPoolTPCnvBase<CaloHitCollection, CaloHitCollection_p1a>
+{
+ public:
+
+ CaloHitCollectionCnv_p1a() {};
+
+ virtual CaloHitCollection* createTransient(const CaloHitCollection_p1a* persObj, MsgStream &log);
+
+ virtual void persToTrans(const CaloHitCollection_p1a* persCont,
+ CaloHitCollection* transCont,
+ MsgStream &log) ;
+ virtual void transToPers(const CaloHitCollection* transCont,
+ CaloHitCollection_p1a* persCont,
+ MsgStream &log) ;
+
+ private:
+
+ static const double m_persEneUnit;
+ static const double m_persLenUnit;
+ static const double m_persAngUnit;
+ static const double m_2bHalfMaximum;
+ static const int m_2bMaximum;
+};
+
+#endif
diff --git a/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1a.h b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1a.h
new file mode 100644
index 000000000..7dc511d9d
--- /dev/null
+++ b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1a.h
@@ -0,0 +1,59 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef CALOHITCOLLECTION_P1A_H
+#define CALOHITCOLLECTION_P1A_H
+
+/*
+
+Authors: Davide Costanzo Rob Duxfield
+Modified for FASER by D. Casper - change internal unsigned short counters to unsigned long
+p2, p3, etc reserved for possible future ATLAS changes
+
+*/
+
+#include <vector>
+#include <string>
+
+class CaloHitCollection_p1a
+{
+ public:
+/// Default constructor
+ CaloHitCollection_p1a ();
+ //private:
+
+ std::vector<float> m_hit1_meanTime; // 1 element per string
+ std::vector<float> m_hit1_x0; //
+ std::vector<float> m_hit1_y0; //
+ std::vector<float> m_hit1_z0; //
+ std::vector<float> m_hit1_theta; //
+ std::vector<float> m_hit1_phi; //
+ std::vector<unsigned long> m_nHits; //
+
+ std::vector<unsigned short> m_hitEne_2b; // 1 element per hit
+ std::vector<unsigned short> m_hitLength_2b; //
+
+ std::vector<unsigned short> m_dTheta; // 1 element per hit except for first hit in string
+ std::vector<unsigned short> m_dPhi; //
+
+ std::vector<float> m_hitEne_4b; // 1 element per hit with m_hitEne_2b[i] == 2**16
+
+ std::vector<float> m_hitLength_4b; // 1 element per hit with m_hitLength_2b[i] == 2**16
+
+ std::vector<unsigned long> m_barcode;
+ std::vector<unsigned long> m_mcEvtIndex;
+ std::vector<char> m_evtColl;
+ std::vector<unsigned long> m_nBC;
+
+ std::vector<unsigned long> m_id;
+ std::vector<unsigned long> m_nId;
+};
+
+
+// inlines
+
+inline
+CaloHitCollection_p1a::CaloHitCollection_p1a () {}
+
+#endif
diff --git a/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloSimEventTPCnvDict.h b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloSimEventTPCnvDict.h
index 80937ffc8..de24c20ac 100644
--- a/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloSimEventTPCnvDict.h
+++ b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/CaloSimEventTPCnvDict.h
@@ -15,6 +15,8 @@
#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCnv_p1.h"
#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1.h"
#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1.h"
+#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1a.h"
+#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1a.h"
#include "FaserCaloSimEventTPCnv/CaloHits/CaloHit_p1.h"
#endif // CALOEVENTTPCNV_CALOSIMEVENTTPCNVDICT_H
diff --git a/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/selection.xml b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/selection.xml
index d73b5ce1a..e90551988 100644
--- a/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/selection.xml
+++ b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/FaserCaloSimEventTPCnv/selection.xml
@@ -4,4 +4,6 @@
<class name="CaloHit_p1" />
<class name="std::vector<CaloHit_p1>" />
<class name="CaloHitCollection_p1" id="134E8045-AB99-43EF-9AD1-324C48830B64" />
+ <class name="CaloHitCollection_p1a" id="02C108EE-0AD9-4444-83BD-D5273FCDEF6F" />
+
</lcgdict>
diff --git a/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/src/CaloHits/CaloHitCollectionCnv_p1a.cxx b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/src/CaloHits/CaloHitCollectionCnv_p1a.cxx
new file mode 100644
index 000000000..d1d6ec53b
--- /dev/null
+++ b/Calorimeter/CaloEventCnv/FaserCaloSimEventTPCnv/src/CaloHits/CaloHitCollectionCnv_p1a.cxx
@@ -0,0 +1,339 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "FaserCaloSimEvent/CaloHit.h"
+#include "FaserCaloSimEvent/CaloHitCollection.h"
+#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollection_p1a.h"
+#include "FaserCaloSimEventTPCnv/CaloHits/CaloHitCollectionCnv_p1a.h"
+#include "GeneratorObjects/HepMcParticleLink.h"
+
+#include <cmath>
+
+//CLHEP
+#include "CLHEP/Geometry/Point3D.h"
+// Gaudi
+#include "GaudiKernel/MsgStream.h"
+#include "GaudiKernel/ThreadLocalContext.h"
+
+// Athena
+#include "AthenaKernel/ExtendedEventContext.h"
+#include "StoreGate/StoreGateSvc.h"
+
+// * * * stolen from eflowRec * * * //
+inline double phicorr(double a)
+{
+ if (a <= -M_PI)
+ {
+ return a+(2*M_PI*floor(-(a-M_PI)/(2*M_PI)));
+ }
+ else if (a > M_PI)
+ {
+ return a-(2*M_PI*floor((a+M_PI)/(2*M_PI)));
+ }
+ else
+ {
+ return a;
+ }
+}
+
+// * * * stolen from eflowRec * * * //
+inline double cycle(double a, double b)
+{
+ double del = b-a;
+ if (del > M_PI)
+ {
+ return a+2.0*M_PI;
+ }
+ else if (del < -M_PI)
+ {
+ return a-2.0*M_PI;
+ }
+ else
+ {
+ return a;
+ }
+}
+
+
+const double CaloHitCollectionCnv_p1a::m_persEneUnit = 1.0e-5;
+const double CaloHitCollectionCnv_p1a::m_persLenUnit = 1.0e-5;
+const double CaloHitCollectionCnv_p1a::m_persAngUnit = 1.0e-5;
+const double CaloHitCollectionCnv_p1a::m_2bHalfMaximum = pow(2.0, 15.0);
+const int CaloHitCollectionCnv_p1a::m_2bMaximum = (unsigned short)(-1);
+
+
+void CaloHitCollectionCnv_p1a::transToPers(const CaloHitCollection* transCont, CaloHitCollection_p1a* persCont, MsgStream &log)
+{
+ // Finds hits belonging to a "string" (in which the end point of one hit is the same as the start point of the next) and
+ // persistifies the end point of each hit plus the start point of the first hit in each string.
+ //
+ // Further compression is achieved by optimising the storage of the position vectors:- start (x,y,z) and (theta,phi) of
+ // first hit are stored as floats, (delta_theta,delta_phi) relative to the fisrst hit are stored as 2 byte numbers and
+ // used to specify the hit direction. All hit lengths are stored as 2 byte numbers.
+ //
+ // Additional savings are achieved by storing the energy loss for each hit as a 2 byte number and only storing the mean
+ // time of the first hit per string.
+ //
+ // See http://indico.cern.ch/getFile.py/access?contribId=11&resId=2&materialId=slides&confId=30893 for more info.
+
+ static const double dRcut = 1.0e-7;
+ static const double dTcut = 1.0;
+
+ const EventContext& ctx = Gaudi::Hive::currentContext();
+ const IProxyDict* proxy = Atlas::getExtendedEventContext(ctx).proxy();
+ const HepMcParticleLink * lastLink=nullptr;
+ int lastId = -1;
+ double stringFirstTheta = 0.0;
+ double stringFirstPhi = 0.0;
+ double lastT = 0.0;
+ double persSumE = 0.0;
+ double transSumE = 0.0;
+ unsigned int idx = 0;
+ unsigned int endBC = 0;
+ unsigned int endId = 0;
+ unsigned int endHit = 0;
+ HepGeom::Point3D<double> lastTransEnd(0.0, 0.0, 0.0);
+ HepGeom::Point3D<double> lastPersEnd(0.0, 0.0, 0.0);
+
+ for (CaloHitCollection::const_iterator it = transCont->begin(); it != transCont->end(); ++it) {
+
+ CaloHitCollection::const_iterator siHit = it;
+
+
+ if ( !lastLink || (siHit->particleLink() != *lastLink) ) {
+
+ // store barcode, eventIndex and McEventCollection once for set of consecutive hits with same barcode
+
+ lastLink = &(siHit->particleLink());
+ persCont->m_barcode.push_back(lastLink->barcode());
+ unsigned short index{0};
+ const HepMcParticleLink::index_type position =
+ HepMcParticleLink::getEventPositionInCollection(lastLink->eventIndex(),
+ lastLink->getEventCollection(),
+ proxy).at(0);
+ if (position!=0) {
+ index = lastLink->eventIndex();
+ if(lastLink->eventIndex()!=static_cast<HepMcParticleLink::index_type>(index)) {
+ log << MSG::WARNING << "Attempting to persistify an eventIndex larger than max unsigned short!" << endmsg;
+ }
+ }
+ persCont->m_mcEvtIndex.push_back(index);
+ persCont->m_evtColl.push_back(lastLink->getEventCollectionAsChar());
+
+ if (idx > 0) {
+ persCont->m_nBC.push_back(idx - endBC);
+ endBC = idx;
+ }
+ }
+
+ if ( (int)siHit->identify() != lastId ) {
+
+ // store id once for set of consecutive hits with same barcode
+
+ lastId = siHit->identify();
+ persCont->m_id.push_back(lastId);
+
+ if (idx > 0) {
+ persCont->m_nId.push_back(idx - endId);
+ endId = idx;
+ }
+ }
+
+ HepGeom::Point3D<double> st = siHit->localStartPosition();
+ HepGeom::Point3D<double> en = siHit->localEndPosition();
+
+ const double dx = st.x() - lastTransEnd.x();
+ const double dy = st.y() - lastTransEnd.y();
+ const double dz = st.z() - lastTransEnd.z();
+ const double t = siHit->meanTime();
+
+ const double dRLast = sqrt(dx * dx + dy * dy + dz * dz); // dR between end of previous hit and start of current one
+ const double dTLast = fabs(t - lastT);
+
+ CLHEP::Hep3Vector direction(0.0, 0.0, 0.0);
+ double theta = 0.0;
+ double phi = 0.0;
+ bool startNewString = false;
+
+ if (dRLast < dRcut && dTLast < dTcut) {
+
+ // hit is part of existing string
+
+ direction = CLHEP::Hep3Vector( en.x() - lastPersEnd.x(), en.y() - lastPersEnd.y(), en.z() - lastPersEnd.z() );
+
+ theta = direction.theta();
+ phi = phicorr( direction.phi() );
+
+ const int dTheta_2b = (int)( (theta - stringFirstTheta) / m_persAngUnit + m_2bHalfMaximum + 0.5 );
+ const int dPhi_2b = (int)( (cycle(phi, stringFirstPhi) - stringFirstPhi) / m_persAngUnit + m_2bHalfMaximum + 0.5 );
+
+ if ( dTheta_2b < m_2bMaximum && dTheta_2b >= 0 && dPhi_2b < m_2bMaximum && dPhi_2b >= 0) {
+ persCont->m_dTheta.push_back(dTheta_2b);
+ persCont->m_dPhi.push_back(dPhi_2b);
+ theta = stringFirstTheta + ( (double)dTheta_2b - m_2bHalfMaximum ) * m_persAngUnit;
+ phi = stringFirstPhi + ( (double)dPhi_2b - m_2bHalfMaximum ) * m_persAngUnit;
+ phi = phicorr(phi);
+ }
+ else {
+ startNewString = true;
+ }
+ }
+
+ if (startNewString || dRLast >= dRcut || dTLast >= dTcut) {
+
+ // begin new hit string
+
+ direction = CLHEP::Hep3Vector( en.x() - st.x(), en.y() - st.y(), en.z() - st.z() );
+
+ theta = direction.theta();
+ phi = phicorr( direction.phi() );
+
+ persCont->m_hit1_meanTime.push_back(t);
+ persCont->m_hit1_x0.push_back(st.x());
+ persCont->m_hit1_y0.push_back(st.y());
+ persCont->m_hit1_z0.push_back(st.z());
+ persCont->m_hit1_theta.push_back(theta);
+ persCont->m_hit1_phi.push_back(phi);
+
+ lastPersEnd = HepGeom::Point3D<double>(st.x(), st.y(), st.z());
+
+ stringFirstTheta = theta;
+ stringFirstPhi = phi;
+
+ if (idx > 0) {
+ persCont->m_nHits.push_back(idx - endHit);
+ endHit = idx;
+ }
+ }
+
+ lastTransEnd = HepGeom::Point3D<double>(en.x(), en.y(), en.z());
+ transSumE += siHit->energyLoss();
+
+ const int eneLoss_2b = (int)((transSumE - persSumE) / m_persEneUnit + 0.5); // calculated to allow recovery sum over
+ // whole hit string to chosen precision
+
+ const int hitLength_2b = (int)(direction.mag() / m_persLenUnit + 0.5); // calculated to give the correct position to
+ // the chosen precision, NOT the length of the
+ // hit (small difference in practice).
+ double eneLoss = 0.0;
+
+ if (eneLoss_2b >= m_2bMaximum) {
+ eneLoss = siHit->energyLoss();
+ persCont->m_hitEne_2b.push_back(m_2bMaximum);
+ persCont->m_hitEne_4b.push_back(eneLoss);
+ }
+ else {
+ eneLoss = eneLoss_2b * m_persEneUnit;
+ persCont->m_hitEne_2b.push_back(eneLoss_2b);
+ }
+
+ double length = 0.0;
+
+ if (hitLength_2b >= m_2bMaximum) {
+ length = direction.mag();
+ persCont->m_hitLength_2b.push_back(m_2bMaximum);
+ persCont->m_hitLength_4b.push_back(direction.mag());
+ }
+ else {
+ length = hitLength_2b * m_persLenUnit;
+ persCont->m_hitLength_2b.push_back(hitLength_2b);
+ }
+
+ CLHEP::Hep3Vector persDir(length, 0.0, 0.0);
+ persDir.setTheta(theta);
+ persDir.setPhi(phi);
+
+ lastPersEnd = (CLHEP::Hep3Vector)lastPersEnd + persDir;
+ persSumE += eneLoss;
+ lastT = t;
+
+ ++idx;
+ }
+
+ persCont->m_nBC.push_back(idx - endBC);
+ persCont->m_nId.push_back(idx - endId);
+ persCont->m_nHits.push_back(idx - endHit);
+}
+
+
+CaloHitCollection* CaloHitCollectionCnv_p1a::createTransient(const CaloHitCollection_p1a* persObj, MsgStream &log) {
+ for (size_t i = 0; i < persObj->m_evtColl.size(); i++)
+ {
+ if (persObj->m_evtColl[i] != 'a')
+ log << MSG::ALWAYS << "Corrupted in createTransient: " << persObj->m_evtColl[i] << endmsg;
+ }
+ std::unique_ptr<CaloHitCollection> trans(std::make_unique<CaloHitCollection>("DefaultCollectionName",persObj->m_nHits.size()));
+ persToTrans(persObj, trans.get(), log);
+ return(trans.release());
+}
+
+
+void CaloHitCollectionCnv_p1a::persToTrans(const CaloHitCollection_p1a* persCont, CaloHitCollection* transCont, MsgStream &/*log*/)
+{
+ const EventContext& ctx = Gaudi::Hive::currentContext();
+
+ unsigned int hitCount = 0;
+ unsigned int angleCount = 0;
+ unsigned int idxBC = 0;
+ unsigned int idxId = 0;
+ unsigned int idxEne4b = 0;
+ unsigned int idxLen4b = 0;
+ unsigned int endHit = 0;
+ unsigned int endBC = 0;
+ unsigned int endId = 0;
+
+
+ for (unsigned int i = 0; i < persCont->m_nHits.size(); i++) {
+
+ if (persCont->m_nHits[i]) {
+ const unsigned int start = endHit;
+ endHit += persCont->m_nHits[i];
+
+ const double t0 = persCont->m_hit1_meanTime[i];
+ const double theta0 = persCont->m_hit1_theta[i];
+ const double phi0 = persCont->m_hit1_phi[i];
+ HepGeom::Point3D<double> endLast(persCont->m_hit1_x0[i], persCont->m_hit1_y0[i], persCont->m_hit1_z0[i]);
+
+ for (unsigned int j = start; j < endHit; j++) {
+ if (j >= endBC + persCont->m_nBC[idxBC])
+ endBC += persCont->m_nBC[idxBC++];
+
+ if (j >= endId + persCont->m_nId[idxId])
+ endId += persCont->m_nId[idxId++];
+
+ const double eneLoss_2b = persCont->m_hitEne_2b[hitCount];
+ const double hitLength_2b = persCont->m_hitLength_2b[hitCount];
+
+ const double eneLoss = (eneLoss_2b < m_2bMaximum) ? eneLoss_2b * m_persEneUnit : persCont->m_hitEne_4b[idxEne4b++];
+ const double length = (hitLength_2b < m_2bMaximum) ? hitLength_2b * m_persLenUnit : persCont->m_hitLength_4b[idxLen4b++];
+
+ const double dTheta = (j > start) ? ((double)persCont->m_dTheta[angleCount] - m_2bHalfMaximum) * m_persAngUnit : 0.0;
+ const double dPhi = (j > start) ? ((double)persCont->m_dPhi[angleCount] - m_2bHalfMaximum) * m_persAngUnit : 0.0;
+
+ const double meanTime = t0;
+ const double theta = theta0 + dTheta;
+ const double phi = phicorr(phi0 + dPhi);
+
+ CLHEP::Hep3Vector r(length, 0.0, 0.0);
+ r.setTheta(theta);
+ r.setPhi(phi);
+
+ HepGeom::Point3D<double> endThis( endLast + r );
+
+ HepMcParticleLink::PositionFlag flag = HepMcParticleLink::IS_INDEX;
+ if (persCont->m_mcEvtIndex[idxBC] == 0) {
+ flag = HepMcParticleLink::IS_POSITION;
+ }
+
+ HepMcParticleLink partLink( persCont->m_barcode[idxBC], persCont->m_mcEvtIndex[idxBC], HepMcParticleLink::ExtendedBarCode::eventCollectionFromChar(persCont->m_evtColl[idxBC]), flag, ctx );
+ transCont->Emplace( endLast, endThis, eneLoss, meanTime, partLink, persCont->m_id[idxId]);
+
+ endLast = endThis;
+
+ ++hitCount;
+ if (j > start) ++angleCount;
+ }
+ }
+ }
+}
diff --git a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.cxx b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.cxx
index 36ff6569e..715a49b37 100644
--- a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.cxx
+++ b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.cxx
@@ -4,6 +4,7 @@
#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1.h"
#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHit_p1.h"
+#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h"
#include "NeutrinoHitCollectionCnv.h"
@@ -18,12 +19,16 @@ NeutrinoHitCollection_PERS* NeutrinoHitCollectionCnv::createPersistent(NeutrinoH
NeutrinoHitCollection* NeutrinoHitCollectionCnv::createTransient() {
MsgStream mlog(msgSvc(), "NeutrinoHitCollectionConverter" );
NeutrinoHitCollectionCnv_p1 converter_p1;
+ NeutrinoHitCollectionCnv_p1a converter_p1a;
static const pool::Guid p1_guid("2CA7AF71-1494-4378-BED4-AFB5C54398AA");
- // static const pool::Guid p1_guid("B2573A16-4B46-4E1E-98E3-F93421680779");
+ static const pool::Guid p1a_guid("0A3CD37D-64CE-405D-98CF-595D678B14B7");
NeutrinoHitCollection *trans_cont(0);
- if( this->compareClassGuid(p1_guid)) {
+ if( this->compareClassGuid(p1a_guid)) {
+ std::unique_ptr< NeutrinoHitCollection_p1a > col_vect( this->poolReadObject< NeutrinoHitCollection_p1a >() );
+ trans_cont = converter_p1a.createTransient( col_vect.get(), mlog );
+ } else if( this->compareClassGuid(p1_guid)) {
std::unique_ptr< NeutrinoHitCollection_p1 > col_vect( this->poolReadObject< NeutrinoHitCollection_p1 >() );
trans_cont = converter_p1.createTransient( col_vect.get(), mlog );
} else {
diff --git a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.h b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.h
index 5ccd3e1bb..a8158c01a 100644
--- a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.h
+++ b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventAthenaPool/src/NeutrinoHitCollectionCnv.h
@@ -6,14 +6,14 @@
#define NEUTRINOHITCOLLECTIONCNV
#include "NeutrinoSimEvent/NeutrinoHitCollection.h"
-#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1.h"
-#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1.h"
+#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1a.h"
+#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h"
#include "AthenaPoolCnvSvc/T_AthenaPoolCustomCnv.h"
// Gaudi
#include "GaudiKernel/MsgStream.h"
// typedef to the latest persistent version
-typedef NeutrinoHitCollection_p1 NeutrinoHitCollection_PERS;
-typedef NeutrinoHitCollectionCnv_p1 NeutrinoHitCollectionCnv_PERS;
+typedef NeutrinoHitCollection_p1a NeutrinoHitCollection_PERS;
+typedef NeutrinoHitCollectionCnv_p1a NeutrinoHitCollectionCnv_PERS;
class NeutrinoHitCollectionCnv : public T_AthenaPoolCustomCnv<NeutrinoHitCollection, NeutrinoHitCollection_PERS > {
friend class CnvFactory<NeutrinoHitCollectionCnv>;
diff --git a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h
new file mode 100644
index 000000000..5666f66e7
--- /dev/null
+++ b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h
@@ -0,0 +1,41 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef NEUTRINOHITCOLLECTIONCNV_P1A_H
+#define NEUTRINOHITCOLLECTIONCNV_P1A_H
+
+// NeutrinoHitCollectionCnv_p1a, T/P separation of Neutrino Hits
+// author D.Costanzo <davide.costanzo@cern.ch>
+// author O.Arnaez <olivier.arnaez@cern.ch>
+
+#include "AthenaPoolCnvSvc/T_AthenaPoolTPConverter.h"
+#include "NeutrinoSimEvent/NeutrinoHitCollection.h"
+#include "NeutrinoHitCollection_p1a.h"
+
+
+class NeutrinoHitCollectionCnv_p1a : public T_AthenaPoolTPCnvBase<NeutrinoHitCollection, NeutrinoHitCollection_p1a>
+{
+ public:
+
+ NeutrinoHitCollectionCnv_p1a() {};
+
+ virtual NeutrinoHitCollection* createTransient(const NeutrinoHitCollection_p1a* persObj, MsgStream &log);
+
+ virtual void persToTrans(const NeutrinoHitCollection_p1a* persCont,
+ NeutrinoHitCollection* transCont,
+ MsgStream &log) ;
+ virtual void transToPers(const NeutrinoHitCollection* transCont,
+ NeutrinoHitCollection_p1a* persCont,
+ MsgStream &log) ;
+
+ private:
+
+ static const double m_persEneUnit;
+ static const double m_persLenUnit;
+ static const double m_persAngUnit;
+ static const double m_2bHalfMaximum;
+ static const int m_2bMaximum;
+};
+
+#endif
diff --git a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1a.h b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1a.h
new file mode 100644
index 000000000..639ff9840
--- /dev/null
+++ b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1a.h
@@ -0,0 +1,57 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef NEUTRINOHITCOLLECTION_P1A_H
+#define NEUTRINOHITCOLLECTION_P1A_H
+
+/*
+
+Authors: Davide Costanzo Rob Duxfield
+
+*/
+
+#include <vector>
+#include <string>
+
+class NeutrinoHitCollection_p1a
+{
+ public:
+/// Default constructor
+ NeutrinoHitCollection_p1a ();
+ //private:
+
+ std::vector<float> m_hit1_meanTime; // 1 element per string
+ std::vector<float> m_hit1_x0; //
+ std::vector<float> m_hit1_y0; //
+ std::vector<float> m_hit1_z0; //
+ std::vector<float> m_hit1_theta; //
+ std::vector<float> m_hit1_phi; //
+ std::vector<unsigned long> m_nHits; //
+
+ std::vector<unsigned short> m_hitEne_2b; // 1 element per hit
+ std::vector<unsigned short> m_hitLength_2b; //
+
+ std::vector<unsigned short> m_dTheta; // 1 element per hit except for first hit in string
+ std::vector<unsigned short> m_dPhi; //
+
+ std::vector<float> m_hitEne_4b; // 1 element per hit with m_hitEne_2b[i] == 2**16
+
+ std::vector<float> m_hitLength_4b; // 1 element per hit with m_hitLength_2b[i] == 2**16
+
+ std::vector<unsigned long> m_barcode;
+ std::vector<unsigned short> m_mcEvtIndex;
+ std::vector<char> m_evtColl;
+ std::vector<unsigned long> m_nBC;
+
+ std::vector<unsigned long> m_id;
+ std::vector<unsigned long> m_nId;
+};
+
+
+// inlines
+
+inline
+NeutrinoHitCollection_p1a::NeutrinoHitCollection_p1a () {}
+
+#endif
diff --git a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnvDict.h b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnvDict.h
index 064ad71ac..15b7afc5b 100644
--- a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnvDict.h
+++ b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnvDict.h
@@ -15,6 +15,9 @@
#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCnv_p1.h"
#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1.h"
#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1.h"
+#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h"
+#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1a.h"
+
#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHit_p1.h"
#endif // NEUTRINOEVENTTPCNV_NEUTRINOSIMEVENTTPCNVDICT_H
diff --git a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/selection.xml b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/selection.xml
index 6808c217a..3332a0a0c 100644
--- a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/selection.xml
+++ b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/NeutrinoSimEventTPCnv/selection.xml
@@ -4,4 +4,5 @@
<class name="NeutrinoHit_p1" />
<class name="std::vector<NeutrinoHit_p1>" />
<class name="NeutrinoHitCollection_p1" id="2CA7AF71-1494-4378-BED4-AFB5C54398AA" />
+ <class name="NeutrinoHitCollection_p1a" id="0A3CD37D-64CE-405D-98CF-595D678B14B7" />
</lcgdict>
diff --git a/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/src/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.cxx b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/src/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.cxx
new file mode 100644
index 000000000..6ef0bb107
--- /dev/null
+++ b/Neutrino/NeutrinoEventCnv/NeutrinoSimEventTPCnv/src/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.cxx
@@ -0,0 +1,334 @@
+/*
+ Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "NeutrinoSimEvent/NeutrinoHit.h"
+#include "NeutrinoSimEvent/NeutrinoHitCollection.h"
+#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollection_p1a.h"
+#include "NeutrinoSimEventTPCnv/NeutrinoHits/NeutrinoHitCollectionCnv_p1a.h"
+#include "GeneratorObjects/HepMcParticleLink.h"
+
+#include <cmath>
+
+//CLHEP
+#include "CLHEP/Geometry/Point3D.h"
+// Gaudi
+#include "GaudiKernel/MsgStream.h"
+#include "GaudiKernel/ThreadLocalContext.h"
+
+// Athena
+#include "AthenaKernel/ExtendedEventContext.h"
+#include "StoreGate/StoreGateSvc.h"
+
+// * * * stolen from eflowRec * * * //
+inline double phicorr(double a)
+{
+ if (a <= -M_PI)
+ {
+ return a+(2*M_PI*floor(-(a-M_PI)/(2*M_PI)));
+ }
+ else if (a > M_PI)
+ {
+ return a-(2*M_PI*floor((a+M_PI)/(2*M_PI)));
+ }
+ else
+ {
+ return a;
+ }
+}
+
+// * * * stolen from eflowRec * * * //
+inline double cycle(double a, double b)
+{
+ double del = b-a;
+ if (del > M_PI)
+ {
+ return a+2.0*M_PI;
+ }
+ else if (del < -M_PI)
+ {
+ return a-2.0*M_PI;
+ }
+ else
+ {
+ return a;
+ }
+}
+
+
+const double NeutrinoHitCollectionCnv_p1a::m_persEneUnit = 1.0e-5;
+const double NeutrinoHitCollectionCnv_p1a::m_persLenUnit = 1.0e-5;
+const double NeutrinoHitCollectionCnv_p1a::m_persAngUnit = 1.0e-5;
+const double NeutrinoHitCollectionCnv_p1a::m_2bHalfMaximum = pow(2.0, 15.0);
+const int NeutrinoHitCollectionCnv_p1a::m_2bMaximum = (unsigned short)(-1);
+
+
+void NeutrinoHitCollectionCnv_p1a::transToPers(const NeutrinoHitCollection* transCont, NeutrinoHitCollection_p1a* persCont, MsgStream &log)
+{
+ // Finds hits belonging to a "string" (in which the end point of one hit is the same as the start point of the next) and
+ // persistifies the end point of each hit plus the start point of the first hit in each string.
+ //
+ // Further compression is achieved by optimising the storage of the position vectors:- start (x,y,z) and (theta,phi) of
+ // first hit are stored as floats, (delta_theta,delta_phi) relative to the fisrst hit are stored as 2 byte numbers and
+ // used to specify the hit direction. All hit lengths are stored as 2 byte numbers.
+ //
+ // Additional savings are achieved by storing the energy loss for each hit as a 2 byte number and only storing the mean
+ // time of the first hit per string.
+ //
+ // See http://indico.cern.ch/getFile.py/access?contribId=11&resId=2&materialId=slides&confId=30893 for more info.
+
+ static const double dRcut = 1.0e-7;
+ static const double dTcut = 1.0;
+
+ const EventContext& ctx = Gaudi::Hive::currentContext();
+ const IProxyDict* proxy = Atlas::getExtendedEventContext(ctx).proxy();
+ const HepMcParticleLink * lastLink=nullptr;
+ int lastId = -1;
+ double stringFirstTheta = 0.0;
+ double stringFirstPhi = 0.0;
+ double lastT = 0.0;
+ double persSumE = 0.0;
+ double transSumE = 0.0;
+ unsigned int idx = 0;
+ unsigned int endBC = 0;
+ unsigned int endId = 0;
+ unsigned int endHit = 0;
+ HepGeom::Point3D<double> lastTransEnd(0.0, 0.0, 0.0);
+ HepGeom::Point3D<double> lastPersEnd(0.0, 0.0, 0.0);
+
+ for (NeutrinoHitCollection::const_iterator it = transCont->begin(); it != transCont->end(); ++it) {
+
+ NeutrinoHitCollection::const_iterator siHit = it;
+
+
+ if ( !lastLink || (siHit->particleLink() != *lastLink) ) {
+
+ // store barcode once for set of consecutive hits with same barcode
+
+ lastLink = &(siHit->particleLink());
+ persCont->m_barcode.push_back(lastLink->barcode());
+ unsigned short index{0};
+ const HepMcParticleLink::index_type position =
+ HepMcParticleLink::getEventPositionInCollection(lastLink->eventIndex(),
+ lastLink->getEventCollection(),
+ proxy).at(0);
+ if (position!=0) {
+ index = lastLink->eventIndex();
+ if(lastLink->eventIndex()!=static_cast<HepMcParticleLink::index_type>(index)) {
+ log << MSG::WARNING << "Attempting to persistify an eventIndex larger than max unsigned short!" << endmsg;
+ }
+ }
+ persCont->m_mcEvtIndex.push_back(index);
+ persCont->m_evtColl.push_back(lastLink->getEventCollectionAsChar());
+
+ if (idx > 0) {
+ persCont->m_nBC.push_back(idx - endBC);
+ endBC = idx;
+ }
+ }
+
+ if ( (int)siHit->identify() != lastId ) {
+
+ // store barcode, eventIndex and McEventCollection once for set of consecutive hits with same barcode
+
+ lastId = siHit->identify();
+ persCont->m_id.push_back(lastId);
+
+ if (idx > 0) {
+ persCont->m_nId.push_back(idx - endId);
+ endId = idx;
+ }
+ }
+
+ HepGeom::Point3D<double> st = siHit->localStartPosition();
+ HepGeom::Point3D<double> en = siHit->localEndPosition();
+
+ const double dx = st.x() - lastTransEnd.x();
+ const double dy = st.y() - lastTransEnd.y();
+ const double dz = st.z() - lastTransEnd.z();
+ const double t = siHit->meanTime();
+
+ const double dRLast = sqrt(dx * dx + dy * dy + dz * dz); // dR between end of previous hit and start of current one
+ const double dTLast = fabs(t - lastT);
+
+ CLHEP::Hep3Vector direction(0.0, 0.0, 0.0);
+ double theta = 0.0;
+ double phi = 0.0;
+ bool startNewString = false;
+
+ if (dRLast < dRcut && dTLast < dTcut) {
+
+ // hit is part of existing string
+
+ direction = CLHEP::Hep3Vector( en.x() - lastPersEnd.x(), en.y() - lastPersEnd.y(), en.z() - lastPersEnd.z() );
+
+ theta = direction.theta();
+ phi = phicorr( direction.phi() );
+
+ const int dTheta_2b = (int)( (theta - stringFirstTheta) / m_persAngUnit + m_2bHalfMaximum + 0.5 );
+ const int dPhi_2b = (int)( (cycle(phi, stringFirstPhi) - stringFirstPhi) / m_persAngUnit + m_2bHalfMaximum + 0.5 );
+
+ if ( dTheta_2b < m_2bMaximum && dTheta_2b >= 0 && dPhi_2b < m_2bMaximum && dPhi_2b >= 0) {
+ persCont->m_dTheta.push_back(dTheta_2b);
+ persCont->m_dPhi.push_back(dPhi_2b);
+ theta = stringFirstTheta + ( (double)dTheta_2b - m_2bHalfMaximum ) * m_persAngUnit;
+ phi = stringFirstPhi + ( (double)dPhi_2b - m_2bHalfMaximum ) * m_persAngUnit;
+ phi = phicorr(phi);
+ }
+ else {
+ startNewString = true;
+ }
+ }
+
+ if (startNewString || dRLast >= dRcut || dTLast >= dTcut) {
+
+ // begin new hit string
+
+ direction = CLHEP::Hep3Vector( en.x() - st.x(), en.y() - st.y(), en.z() - st.z() );
+
+ theta = direction.theta();
+ phi = phicorr( direction.phi() );
+
+ persCont->m_hit1_meanTime.push_back(t);
+ persCont->m_hit1_x0.push_back(st.x());
+ persCont->m_hit1_y0.push_back(st.y());
+ persCont->m_hit1_z0.push_back(st.z());
+ persCont->m_hit1_theta.push_back(theta);
+ persCont->m_hit1_phi.push_back(phi);
+
+ lastPersEnd = HepGeom::Point3D<double>(st.x(), st.y(), st.z());
+
+ stringFirstTheta = theta;
+ stringFirstPhi = phi;
+
+ if (idx > 0) {
+ persCont->m_nHits.push_back(idx - endHit);
+ endHit = idx;
+ }
+ }
+
+ lastTransEnd = HepGeom::Point3D<double>(en.x(), en.y(), en.z());
+ transSumE += siHit->energyLoss();
+
+ const int eneLoss_2b = (int)((transSumE - persSumE) / m_persEneUnit + 0.5); // calculated to allow recovery sum over
+ // whole hit string to chosen precision
+
+ const int hitLength_2b = (int)(direction.mag() / m_persLenUnit + 0.5); // calculated to give the correct position to
+ // the chosen precision, NOT the length of the
+ // hit (small difference in practice).
+ double eneLoss = 0.0;
+
+ if (eneLoss_2b >= m_2bMaximum) {
+ eneLoss = siHit->energyLoss();
+ persCont->m_hitEne_2b.push_back(m_2bMaximum);
+ persCont->m_hitEne_4b.push_back(eneLoss);
+ }
+ else {
+ eneLoss = eneLoss_2b * m_persEneUnit;
+ persCont->m_hitEne_2b.push_back(eneLoss_2b);
+ }
+
+ double length = 0.0;
+
+ if (hitLength_2b >= m_2bMaximum) {
+ length = direction.mag();
+ persCont->m_hitLength_2b.push_back(m_2bMaximum);
+ persCont->m_hitLength_4b.push_back(direction.mag());
+ }
+ else {
+ length = hitLength_2b * m_persLenUnit;
+ persCont->m_hitLength_2b.push_back(hitLength_2b);
+ }
+
+ CLHEP::Hep3Vector persDir(length, 0.0, 0.0);
+ persDir.setTheta(theta);
+ persDir.setPhi(phi);
+
+ lastPersEnd = (CLHEP::Hep3Vector)lastPersEnd + persDir;
+ persSumE += eneLoss;
+ lastT = t;
+
+ ++idx;
+ }
+
+ persCont->m_nBC.push_back(idx - endBC);
+ persCont->m_nId.push_back(idx - endId);
+ persCont->m_nHits.push_back(idx - endHit);
+}
+
+
+NeutrinoHitCollection* NeutrinoHitCollectionCnv_p1a::createTransient(const NeutrinoHitCollection_p1a* persObj, MsgStream &log) {
+ std::unique_ptr<NeutrinoHitCollection> trans(std::make_unique<NeutrinoHitCollection>("DefaultCollectionName",persObj->m_nHits.size()));
+ persToTrans(persObj, trans.get(), log);
+ return(trans.release());
+}
+
+
+void NeutrinoHitCollectionCnv_p1a::persToTrans(const NeutrinoHitCollection_p1a* persCont, NeutrinoHitCollection* transCont, MsgStream &/*log*/)
+{
+ const EventContext& ctx = Gaudi::Hive::currentContext();
+
+ unsigned int hitCount = 0;
+ unsigned int angleCount = 0;
+ unsigned int idxBC = 0;
+ unsigned int idxId = 0;
+ unsigned int idxEne4b = 0;
+ unsigned int idxLen4b = 0;
+ unsigned int endHit = 0;
+ unsigned int endBC = 0;
+ unsigned int endId = 0;
+
+ for (unsigned int i = 0; i < persCont->m_nHits.size(); i++) {
+
+ if (persCont->m_nHits[i]) {
+
+ const unsigned int start = endHit;
+ endHit += persCont->m_nHits[i];
+
+ const double t0 = persCont->m_hit1_meanTime[i];
+ const double theta0 = persCont->m_hit1_theta[i];
+ const double phi0 = persCont->m_hit1_phi[i];
+ HepGeom::Point3D<double> endLast(persCont->m_hit1_x0[i], persCont->m_hit1_y0[i], persCont->m_hit1_z0[i]);
+
+ for (unsigned int j = start; j < endHit; j++) {
+
+ if (j >= endBC + persCont->m_nBC[idxBC])
+ endBC += persCont->m_nBC[idxBC++];
+
+ if (j >= endId + persCont->m_nId[idxId])
+ endId += persCont->m_nId[idxId++];
+
+ const double eneLoss_2b = persCont->m_hitEne_2b[hitCount];
+ const double hitLength_2b = persCont->m_hitLength_2b[hitCount];
+
+ const double eneLoss = (eneLoss_2b < m_2bMaximum) ? eneLoss_2b * m_persEneUnit : persCont->m_hitEne_4b[idxEne4b++];
+ const double length = (hitLength_2b < m_2bMaximum) ? hitLength_2b * m_persLenUnit : persCont->m_hitLength_4b[idxLen4b++];
+
+ const double dTheta = (j > start) ? ((double)persCont->m_dTheta[angleCount] - m_2bHalfMaximum) * m_persAngUnit : 0.0;
+ const double dPhi = (j > start) ? ((double)persCont->m_dPhi[angleCount] - m_2bHalfMaximum) * m_persAngUnit : 0.0;
+
+ const double meanTime = t0;
+ const double theta = theta0 + dTheta;
+ const double phi = phicorr(phi0 + dPhi);
+
+ CLHEP::Hep3Vector r(length, 0.0, 0.0);
+ r.setTheta(theta);
+ r.setPhi(phi);
+
+ HepGeom::Point3D<double> endThis( endLast + r );
+
+ HepMcParticleLink::PositionFlag flag = HepMcParticleLink::IS_INDEX;
+ if (persCont->m_mcEvtIndex[idxBC] == 0) {
+ flag = HepMcParticleLink::IS_POSITION;
+ }
+ HepMcParticleLink partLink( persCont->m_barcode[idxBC], persCont->m_mcEvtIndex[idxBC], HepMcParticleLink::ExtendedBarCode::eventCollectionFromChar(persCont->m_evtColl[idxBC]), flag, ctx );
+ transCont->Emplace( endLast, endThis, eneLoss, meanTime, partLink, persCont->m_id[idxId]);
+
+ endLast = endThis;
+
+ ++hitCount;
+ if (j > start) ++angleCount;
+ }
+ }
+ }
+}
--
GitLab