Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add TrkStrawHitMC provenance and guard against empty StrawDigiMCs #1291

Merged
merged 13 commits into from
Jul 1, 2024
Merged
Show file tree
Hide file tree
Changes from 12 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
192 changes: 103 additions & 89 deletions CommonMC/src/SelectRecoMC_module.cc
Original file line number Diff line number Diff line change
Expand Up @@ -53,6 +53,7 @@
#include <iostream>
#include <set>
#include <string>
#include <algorithm>

namespace mu2e {
class SelectRecoMC : public art::EDProducer {
Expand Down Expand Up @@ -233,14 +234,17 @@ namespace mu2e {
// find the referenced sim particle
int spref(-1);
auto const& sdmc = sdmcc.at(hit.index()); // bounds-check for security;
for(size_t isp=0;isp < spcc.size(); isp++){
auto const& spc = spcc[isp];
if(sdmc.earlyStrawGasStep()->simParticle() == spc._spp){
spref = isp;
break;
// if mc info is not meaningful, do not try to inspect any SimParticles
if (sdmc.provenance().ContainsSimulation()){
for(size_t isp=0;isp < spcc.size(); isp++){
auto const& spc = spcc[isp];
if(sdmc.earlyStrawGasStep()->simParticle() == spc._spp){
spref = isp;
break;
}
}
if(spref < 0)throw cet::exception("Reco")<<"mu2e::SelectRecoMC: missing index"<< std::endl;
}
if(spref < 0)throw cet::exception("Reco")<<"mu2e::SelectRecoMC: missing index"<< std::endl;
TrkStrawHitMC tshmc;
fillTSHMC(tshmc,hit.index(),spref,sdmc,tracker,srep);
mcseed._tshmcs.push_back(tshmc);
Expand All @@ -250,26 +254,30 @@ namespace mu2e {
void SelectRecoMC::fillUnusedTSHMC( SPCC const& spcc, StrawDigiMCCollection const& sdmcc,
Tracker const& tracker, std::shared_ptr<const StrawResponse>const& srep,
KalSeedMC& mcseed) {
// either keep hits only from the primary or from all contributing
size_t ispmax = _saveallunused? spcc.size() : 1;
// either keep hits only from the primary or all contributing particles,
// if any exist. if they do not, then there is no primary to default to.
size_t ispmax = _saveallunused? spcc.size() : std::min<int>(spcc.size(), 1);
for(size_t isp=0; isp < ispmax; ++isp){
auto const& spc = spcc[isp];
for (size_t isdmc=0; isdmc < sdmcc.size(); isdmc++){
auto const& sdmc = sdmcc[isdmc];
auto const& sgs = *(sdmc.earlyStrawGasStep());
if(sgs.simParticle() == spc._spp){
// search to see if the associated digi is already on the track
bool used(false);
for(auto const& tshmc : mcseed._tshmcs ) {
if(isdmc == tshmc.strawDigiMCIndex()){
used = true;
break;
// if this contains no MC information, then we cannot inspect it further
if (sdmc.provenance().ContainsSimulation()){
auto const& sgs = *(sdmc.earlyStrawGasStep());
if(sgs.simParticle() == spc._spp){
// search to see if the associated digi is already on the track
bool used(false);
for(auto const& tshmc : mcseed._tshmcs ) {
if(isdmc == tshmc.strawDigiMCIndex()){
used = true;
break;
}
}
if(!used){
TrkStrawHitMC tshmc;
fillTSHMC(tshmc,isdmc,isp,sdmc,tracker,srep);
mcseed._tshmcs.push_back(tshmc);
}
}
if(!used){
TrkStrawHitMC tshmc;
fillTSHMC(tshmc,isdmc,isp,sdmc,tracker,srep);
mcseed._tshmcs.push_back(tshmc);
}
}
}
Expand All @@ -278,78 +286,84 @@ namespace mu2e {

void SelectRecoMC::fillTSHMC(TrkStrawHitMC& tshmc, size_t isdmc, size_t isp,StrawDigiMC const& sdmc,
Tracker const& tracker, std::shared_ptr<const StrawResponse>const& srep ) {
// always record the index, to avoid downstream corruptions
tshmc._sdmcindex = isdmc;
tshmc._spindex = isp;
tshmc._energySum = sdmc.triggerEnergySum(sdmc.earlyEnd());
const auto& sgs = *(sdmc.earlyStrawGasStep());
tshmc._cpos = XYZVectorF(sdmc.clusterPosition(sdmc.earlyEnd()));
tshmc._mom = sgs.momentum();
tshmc._time = sgs.time();
if (_onSpill){
tshmc._time = fmod(tshmc._time,_mbtime);
// fix for DAQ wrapping
if(tshmc._time < -_pbtimemc)tshmc._time += _mbtime;
}
tshmc._strawId = sdmc.strawId();
tshmc._earlyend = sdmc.earlyEnd();
// compute the signal propagation time and drift time
double vprop = 2.0*srep->halfPropV(sdmc.strawId(),1000.0*tshmc._energySum);
auto const& straw = tracker.straw(sdmc.strawId());
auto wdir = XYZVectorF(straw.wireDirection());
auto tdir = sgs.momentum().Unit();
double pdist = (straw.wireEnd(sdmc.earlyEnd())-sdmc.clusterPosition(sdmc.earlyEnd())).dot(straw.wireDirection());
tshmc._tprop = fabs(pdist)/vprop;
tshmc._tdrift = sdmc.wireEndTime(sdmc.earlyEnd()) -tshmc._time - tshmc._tprop - _pbtimemc - 2.4; // temporary kludge offset FIXME!
if (_onSpill)
tshmc._tdrift = fmod(tshmc._tdrift,_mbtime);
const static XYZVectorF bdir(0.0,0.0,1.0);
// propagate interpretability of StrawDigiMC to TrkStrawHitMC
tshmc._provenance = sdmc.provenance();
// if there is no MC information at all, leave TrkStrawHitMC empty
if (tshmc._provenance.ContainsSimulation()){
tshmc._spindex = isp;
tshmc._energySum = sdmc.triggerEnergySum(sdmc.earlyEnd());
const auto& sgs = *(sdmc.earlyStrawGasStep());
tshmc._cpos = XYZVectorF(sdmc.clusterPosition(sdmc.earlyEnd()));
tshmc._mom = sgs.momentum();
tshmc._time = sgs.time();
if (_onSpill){
tshmc._time = fmod(tshmc._time,_mbtime);
// fix for DAQ wrapping
if(tshmc._time < -_pbtimemc)tshmc._time += _mbtime;
}
tshmc._strawId = sdmc.strawId();
tshmc._earlyend = sdmc.earlyEnd();
// compute the signal propagation time and drift time
double vprop = 2.0*srep->halfPropV(sdmc.strawId(),1000.0*tshmc._energySum);
auto const& straw = tracker.straw(sdmc.strawId());
auto wdir = XYZVectorF(straw.wireDirection());
auto tdir = sgs.momentum().Unit();
double pdist = (straw.wireEnd(sdmc.earlyEnd())-sdmc.clusterPosition(sdmc.earlyEnd())).dot(straw.wireDirection());
tshmc._tprop = fabs(pdist)/vprop;
tshmc._tdrift = sdmc.wireEndTime(sdmc.earlyEnd()) -tshmc._time - tshmc._tprop - _pbtimemc - 2.4; // temporary kludge offset FIXME!
if (_onSpill)
tshmc._tdrift = fmod(tshmc._tdrift,_mbtime);
const static XYZVectorF bdir(0.0,0.0,1.0);

TwoLinePCA_XYZ wirepca( XYZVectorF(straw.wirePosition()), XYZVectorF(straw.wireDirection()),
tshmc._cpos, tdir );
TwoLinePCA_XYZ strawpca( XYZVectorF(straw.strawPosition()), XYZVectorF(straw.strawDirection()),
tshmc._cpos, tdir );
TwoLinePCA_XYZ wirepca( XYZVectorF(straw.wirePosition()), XYZVectorF(straw.wireDirection()),
tshmc._cpos, tdir );
TwoLinePCA_XYZ strawpca( XYZVectorF(straw.strawPosition()), XYZVectorF(straw.strawDirection()),
tshmc._cpos, tdir );

// vector from middle of wire to POCA
auto wireMid_to_POCA = wirepca.point2() - XYZVectorF(straw.wirePosition());
// longitudinal distance along wire to POCA
tshmc._wireLen = wireMid_to_POCA.Dot(wdir);
// perpedicular vector from wire to POCA
auto POCA_delta = wirepca.point2() - wirepca.point1();
// define DOCA to be negative in direction of track cross wire
auto doca_sign_dir = tdir.Cross(wdir).Unit();
tshmc._wireDOCA = -1*doca_sign_dir.Dot(POCA_delta);
// phi is defined by vector to POCA, with 0 in 'Bdir'
// and +pi/2 in V dir (radially out)
auto Vdir = wdir.Cross(bdir);
// ensure Vdir is pointing radially out
if (Vdir.Dot(XYZVectorF(straw.wirePosition())) < 0.0) Vdir *= -1.0;
tshmc._wirePhi = atan2(POCA_delta.Dot(Vdir),POCA_delta.Dot(bdir));
// rdrift is the expected T2D given tdrift and phi
tshmc._rdrift = srep->strawDrift().T2D(tshmc._tdrift,tshmc._wirePhi);
// wireDot is cos angle between track and wire
tshmc._wireDot = tdir.Dot(wdir);
// wireTau is delta distance perpedicular to the wire from POCA to trigger cluster
// first define unit vector perpendicular to particle track and wire direction
auto track_cross_wire = tdir.Cross(wdir).Unit();
// then this cross wdir is the vector perpendicular to the wire and the vector to POCA
auto wire_cross_delta = wdir.Cross(track_cross_wire);
// relative position of the trigger cluster
auto wireMid_to_cluster = tshmc._cpos-XYZVectorF(straw.wirePosition());
// wireTau is then the relative position of the trigger cluster in this direction
tshmc._wireTau = wireMid_to_cluster.Dot(wire_cross_delta);
// vector from middle of wire to POCA
auto wireMid_to_POCA = wirepca.point2() - XYZVectorF(straw.wirePosition());
// longitudinal distance along wire to POCA
tshmc._wireLen = wireMid_to_POCA.Dot(wdir);
// perpedicular vector from wire to POCA
auto POCA_delta = wirepca.point2() - wirepca.point1();
// define DOCA to be negative in direction of track cross wire
auto doca_sign_dir = tdir.Cross(wdir).Unit();
tshmc._wireDOCA = -1*doca_sign_dir.Dot(POCA_delta);
// phi is defined by vector to POCA, with 0 in 'Bdir'
// and +pi/2 in V dir (radially out)
auto Vdir = wdir.Cross(bdir);
// ensure Vdir is pointing radially out
if (Vdir.Dot(XYZVectorF(straw.wirePosition())) < 0.0) Vdir *= -1.0;
tshmc._wirePhi = atan2(POCA_delta.Dot(Vdir),POCA_delta.Dot(bdir));
// rdrift is the expected T2D given tdrift and phi
tshmc._rdrift = srep->strawDrift().T2D(tshmc._tdrift,tshmc._wirePhi);
// wireDot is cos angle between track and wire
tshmc._wireDot = tdir.Dot(wdir);
// wireTau is delta distance perpedicular to the wire from POCA to trigger cluster
// first define unit vector perpendicular to particle track and wire direction
auto track_cross_wire = tdir.Cross(wdir).Unit();
// then this cross wdir is the vector perpendicular to the wire and the vector to POCA
auto wire_cross_delta = wdir.Cross(track_cross_wire);
// relative position of the trigger cluster
auto wireMid_to_cluster = tshmc._cpos-XYZVectorF(straw.wirePosition());
// wireTau is then the relative position of the trigger cluster in this direction
tshmc._wireTau = wireMid_to_cluster.Dot(wire_cross_delta);

auto sdir = XYZVectorF(straw.strawDirection());
// perpedicular vector from straw center to POCA
auto sPOCA_delta = strawpca.point2() - strawpca.point1();
// define DOCA to be negative in direction of track cross straw
auto sdoca_sign_dir = tdir.Cross(sdir).Unit();
tshmc._strawDOCA = -1*sdoca_sign_dir.Dot(sPOCA_delta);
// phi is defined by vector to POCA, with 0 in 'Bdir'
// and +pi/2 in V dir (radially out)
auto sVdir = sdir.Cross(bdir);
// ensure Vdir is pointing radially out
if (sVdir.Dot(XYZVectorF(straw.strawPosition())) < 0.0) sVdir *= -1.0;
tshmc._strawPhi = atan2(sPOCA_delta.Dot(sVdir),sPOCA_delta.Dot(bdir));
auto sdir = XYZVectorF(straw.strawDirection());
// perpedicular vector from straw center to POCA
auto sPOCA_delta = strawpca.point2() - strawpca.point1();
// define DOCA to be negative in direction of track cross straw
auto sdoca_sign_dir = tdir.Cross(sdir).Unit();
tshmc._strawDOCA = -1*sdoca_sign_dir.Dot(sPOCA_delta);
// phi is defined by vector to POCA, with 0 in 'Bdir'
// and +pi/2 in V dir (radially out)
auto sVdir = sdir.Cross(bdir);
// ensure Vdir is pointing radially out
if (sVdir.Dot(XYZVectorF(straw.strawPosition())) < 0.0) sVdir *= -1.0;
tshmc._strawPhi = atan2(sPOCA_delta.Dot(sVdir),sPOCA_delta.Dot(bdir));
}
}

void SelectRecoMC::fillSDMCI(KalSeedMC const& mcseed, SDIS& sdindices) {
Expand Down
33 changes: 18 additions & 15 deletions Filters/src/StrawDigiMCFilter_module.cc
Original file line number Diff line number Diff line change
Expand Up @@ -66,21 +66,24 @@ namespace mu2e {
// keep count of digis produced by specific particle
std::map<art::Ptr<SimParticle>,unsigned> pmap;
for(auto const& mcdigi : *mcdigis) {
// look at the early end
StrawEnd fend = mcdigi.earlyEnd();
auto const& step = mcdigi.strawGasStep(fend);
art::Ptr<SimParticle> const& sp = step->simParticle();
auto const& mom = step->momentum(); // cast to 3-vector
if(debug_ > 0)std::cout <<"SimParticle PDG = " << sp->pdgId() << " Mom = " << sqrt(mom.mag2()) << std::endl;
bool goodpdg(true);
if(pdgs_.size() > 0)
goodpdg = std::find(pdgs_.begin(),pdgs_.end(),sp->pdgId()) != pdgs_.end();
if(goodpdg && sqrt(mom.mag2()) > minpmom_ && sqrt(mom.mag2()) < maxpmom_ ){
auto mapfnd = pmap.find(sp);
if(mapfnd == pmap.end())
pmap[sp] = 1;
else
++mapfnd->second;
// do not inspect truth info for digis not produced in simulation
if (mcdigi.provenance().ContainsSimulation()){
// look at the early end
StrawEnd fend = mcdigi.earlyEnd();
auto const& step = mcdigi.strawGasStep(fend);
art::Ptr<SimParticle> const& sp = step->simParticle();
auto const& mom = step->momentum(); // cast to 3-vector
if(debug_ > 0)std::cout <<"SimParticle PDG = " << sp->pdgId() << " Mom = " << sqrt(mom.mag2()) << std::endl;
bool goodpdg(true);
if(pdgs_.size() > 0)
goodpdg = std::find(pdgs_.begin(),pdgs_.end(),sp->pdgId()) != pdgs_.end();
if(goodpdg && sqrt(mom.mag2()) > minpmom_ && sqrt(mom.mag2()) < maxpmom_ ){
auto mapfnd = pmap.find(sp);
if(mapfnd == pmap.end())
pmap[sp] = 1;
else
++mapfnd->second;
}
}
}
// check if any single particle generated enough digis. Save All the particles
Expand Down
32 changes: 32 additions & 0 deletions MCDataProducts/inc/DigiProvenance.hh
Original file line number Diff line number Diff line change
@@ -0,0 +1,32 @@
// Flag to designate a digi as having been produced from simulation,
// read from preexisting data, or a hybridization fo the two
// Ed Callaghan, 2024

#ifndef MCDataProducts_DigiProvenance_hh
#define MCDataProducts_DigiProvenance_hh

// stl
#include <map>
#include <string>

// mu2e
#include "Offline/GeneralUtilities/inc/EnumToStringSparse.hh"

namespace mu2e{
// enum equipped with std::string descriptions
class DigiProvenanceDetail{
public:
enum enum_type {unknown=0, Simulation, Mixed, External};
static std::string const& typeName();
static std::map<enum_type, std::string> const& names();
};
using StringedDigiProvenance = EnumToStringSparse<DigiProvenanceDetail>;

class DigiProvenance: public StringedDigiProvenance{
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Three things. First, EnumToStringSparse does not have a virtual d'tor so we should not inherit from it. That's easy to fix, we can give it a virtual d'tor and also add the other 4 rule of 5 functions as =default. ( My style is if the compiler will do the right thing for all rule of 5 functions I just leave them out and add a one line comment to that effect - I see that EnumToStringSparse is so old that the rule of 5 was then the rule of 3.) Go ahead and edit EnumToStringSparse.

The second issue is the CSAID C++ experts strongly advise against using inheritance when containment will do. We do violate this in some places. I think we are OK here.

Third, I think that if you add
typedef StringedDigiProvenance::enum_type enum_type;

(or try the equivalent using would be better ).

Then downstream code can avoid needing the DigiProvenanceDetail:: and you may also be able to avoid the static casts. I am not 100% sure that I have all details right here but I am pretty sure I am close.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

After sleeping on it, there is an different solution that avoids introducing the third type and avoids issue with inheritance, just use a free function:

using DigiProvenance = EnumToStringSparse<DigiProvenanceDetail>;

bool containsSimulation(  DigiProvenance::enum_type  id ){
   return (  id == DigiProvenance::Simulation || id == DigiProvenance::Mixed );
}

Copy link
Contributor Author

@edcallaghan edcallaghan Jun 28, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The free-function approach appears to be the path of least resistance --- it is preferable, IMO, to the extra baggage needed to make another distinct class pass muster. (FWIW, I went with inheritance over containment here because otherwise there are extra implementations to add to the list (i.e. forwarding all of the member functions of EnumToStringSparse<>)).

That being said, I'm seeing issues with serialization of that simpler implementation (using DigiProvenance = EnumToStringSparse<DigiProvenanceDetail>) manifesting as data corruption after deserializing, in which the actual enum_type not being set correctly. Since this should be simpler the scenario, I'm not really sure what's going wrong. I'll have another try at figuring it out, but if I can't then I may need to just push updates to EnumToStringSparse<> to make it "inheritable."

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry for the slow response.

Can you point me to your working code that is failing and tell me what the error message is?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

No worries. I have the changes in a branch at https://github.com/edcallaghan/mu2e-Offline/tree/update_TrkStrawHitMC_provenance_free_functions. If you compile that (should work at head) and run the contents of ~ejc3/public/DigiProvenance, you can reproduce for yourself. (This is a pair of jobs exemplified in .../run.sh, in which we first physics-mix'' [pm] a few NoPrimary [np] events, and then digi-mix'' [dm] a few CEs into them [ce]). You can see the error in .../dm-ce.log, the relevant bit pasted here:

---- EventProcessorFailure BEGIN
  EventProcessor: an exception occurred during current event processing
  ---- ScheduleExecutionFailure BEGIN
    Path: ProcessingStopped.
    ---- ProductNotFound BEGIN
      A request to resolve an Ptr to a product containing items of type: mu2e::StrawGasStep with ProductID 2189702261
      cannot be satisfied because the product cannot be found.
      The productGetter was not set -- are you trying to dereference a Ptr during mixing?
      The above exception was thrown while processing module StrawDigiMCFilter/Triggerable run: 1202 subRun: 0 event: 1
    ---- ProductNotFound END
    Exception going through path TriggerablePath
  ---- ScheduleExecutionFailure END
---- EventProcessorFailure END

But: my statement that the problem is in the serialization of the provenance was confused and, I think, incorrect. After another round of less-sloppy testing, the problem presents with the current version in this PR as well.

I suspect that the actual problem is with the construction of a new StrawDigiMC at https://github.com/Mu2e/Offline/blob/main/TrackerMC/src/StrawDigiBundleCollection.cc#L223, ultimately stemming from the same "scope constraint" on use of art data products as I mention in #1290. That has nothing to do with the implementation of DigiProvenance, so I can push the simpler implementation for that here. If I'm right, then another PR will need to be opened which addresses the invalid Ptrs.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am looking at it now.

Just a wild guess and it's from memory. If your are making two data products, ACollection and BCollection, both in the same module, and if one contains Ptrs into the other, you cannot dereference the Ptrs until the module has returned. So they only work in downstream modules. The reason is that the pointee data product does not exist until after the return ( the move and the return are like a database two phase commit ).

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please make dm-ce.fcl readable by me.

ls -l ~ejc3/public/DigiProvenance/dm-ce.fcl
-rw-------+ 1 ejc3 nobody 4177 Jun 30 15:46 /nashome/e/ejc3/public/DigiProvenance/dm-ce.fcl

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It looks like my wild guess is irrelevant The changes to Filters/src/StrawDigiMCFilter_module.cc look simple enough. So I agree it's likely that the bug is earlier in the workflow.

That's for making the change to the free function and fixing the typo. My reason for suggesting the free function was the same as in your reply.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry, ~ejc3/public/DigiProvenance/dm-ce.fcl is now public-readable.

But, that being said, I think the issue here is actually much dumber than what I was getting at. When "digi mixing," I am mixing in the StrawDigiMCs, each of which contains an art::Ptr<SimParticle> into a SimParticleCollection. But I don't actually mix the SimParticleCollection in, so of course the pointer can't be dereferenced.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Are all of the fixes in now?

public:
DigiProvenance(StringedDigiProvenance);
bool ContainsSimulation() const;
};
}

#endif
5 changes: 5 additions & 0 deletions MCDataProducts/inc/KalSeedMC.hh
Original file line number Diff line number Diff line change
Expand Up @@ -10,12 +10,14 @@
#include "Offline/DataProducts/inc/GenVector.hh"
#include "Offline/DataProducts/inc/StrawEnd.hh"
#include "Offline/DataProducts/inc/VirtualDetectorId.hh"
#include <Offline/GeneralUtilities/inc/EnumToStringSparse.hh>
#include "Offline/RecoDataProducts/inc/KalSeed.hh"
#include "Offline/MCDataProducts/inc/SimParticle.hh"
#include "Offline/MCDataProducts/inc/ProcessCode.hh"
#include "Offline/MCDataProducts/inc/GenId.hh"
#include "Offline/MCDataProducts/inc/MCRelationship.hh"
#include "Offline/MCDataProducts/inc/CaloClusterMC.hh"
#include "Offline/MCDataProducts/inc/DigiProvenance.hh"
#include "art/Framework/Principal/Handle.h"
#include "cetlib/map_vector.h"
#include <Rtypes.h>
Expand Down Expand Up @@ -70,6 +72,7 @@ namespace mu2e {
//
// MC information for TrackStrawHits on this fit
struct TrkStrawHitMC {
TrkStrawHitMC(): _provenance(static_cast<StringedDigiProvenance>(DigiProvenanceDetail::Simulation)) {}
StrawHitIndex strawDigiMCIndex() const { return _sdmcindex; }
StrawHitIndex simPartStubIndex() const { return _spindex; }
StrawId const& strawid() const { return _strawId; }
Expand All @@ -95,6 +98,7 @@ namespace mu2e {
float _wireTau; // threshold cluster distance to the wire along the perpedicular particle path
float _strawDOCA; // signed doca to straw
float _strawPhi; // cylindrical phi from -pi to pi with 0 in Z direction
DigiProvenance _provenance; // origin/validity of MC info object
};

struct KalSeedMC {
Expand All @@ -104,6 +108,7 @@ namespace mu2e {
std::vector<VDStep> const& vdSteps() const { return _vdsteps; }
TrkStrawHitMC const& trkStrawHitMC(size_t index) const { return _tshmcs.at(index); }
SimPartStub const& simParticle(TrkStrawHitMC const& tshmc) const { return simParticle(tshmc.simPartStubIndex()); }
bool ContainsSimulation() const;
// data products
std::vector<SimPartStub> _simps; // associated sim particles, and their relationship
std::vector<TrkStrawHitMC> _tshmcs; // MC info for each TrkStrawHitSeed
Expand Down
21 changes: 5 additions & 16 deletions MCDataProducts/inc/StrawDigiMC.hh
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,7 @@
#include "Offline/DataProducts/inc/StrawId.hh"
#include "Offline/DataProducts/inc/StrawEnd.hh"
#include "Offline/DataProducts/inc/GenVector.hh"
#include "Offline/MCDataProducts/inc/DigiProvenance.hh"
#include "Offline/MCDataProducts/inc/StrawGasStep.hh"
#include "Offline/MCDataProducts/inc/StepPointMC.hh"
// CLHEP includes
Expand All @@ -22,19 +23,7 @@
#include <array>
#include <string>

// Mu2e includes
#include <Offline/GeneralUtilities/inc/EnumToStringSparse.hh>

namespace mu2e {
// enum equipped with std::string descriptions
class StrawDigiProvenanceDetail{
public:
enum enum_type {unknown=0, Simulation, Mixed, External};
static std::string const& typeName();
static std::map<enum_type, std::string> const& names();
};
using StrawDigiProvenance = EnumToStringSparse<StrawDigiProvenanceDetail>;

class StrawDigiMC{

public:
Expand All @@ -46,16 +35,16 @@ namespace mu2e {

StrawDigiMC();
// construct from hitlets
StrawDigiMC(StrawId sid, PA cpos, FA ctime, FA wetime, SGSPA sgs, StrawDigiProvenance::enum_type=StrawDigiProvenance::Simulation);
StrawDigiMC(StrawId sid, PA cpos, FA ctime, FA wetime, SGSPA sgs, DigiProvenance::enum_type=DigiProvenance::Simulation);

// use compuater copy construcors
StrawDigiMC(const StrawDigiMC& rhs, SGSPA sgsp ); // update the Ptrs
StrawDigiMC(const StrawDigiMC& rhs, StrawDigiProvenance::enum_type provenance ); // update validity
StrawDigiMC(const StrawDigiMC& rhs, DigiProvenance::enum_type provenance ); // update validity

// Accessors
StrawId strawId() const { return _strawid; }

StrawDigiProvenance provenance() const { return _provenance; }
DigiProvenance provenance() const { return _provenance; }

SGSP const& strawGasStep(StrawEnd strawend) const { return _sgspa[strawend]; }
SGSPA const& strawGasSteps() const { return _sgspa; }
Expand Down Expand Up @@ -83,7 +72,7 @@ namespace mu2e {
FA _ctime; // times of the trigger clusters
FA _wtime; // times at the wire ends of the signals which fired the TDC.
SGSPA _sgspa; // StrawGasStep that triggered each end
StrawDigiProvenance _provenance; // level of association with any true MC events
DigiProvenance _provenance; // level of association with any MC truth info
};

inline std::ostream& operator<<( std::ostream& ost,
Expand Down
Loading