LArPandoraInput class. More...

#include <LArPandoraInput.h>

Classes
class	Settings
	Settings class. More...

Static Public Member Functions
static void	CreatePandoraHits2D (const art::Event &evt, const Settings &settings, const LArDriftVolumeMap &driftVolumeMap, const HitVector &hitVector, IdToHitMap &idToHitMap)
	Create the Pandora 2D hits from the ART hits. More...

static void	CreatePandoraLArTPCs (const Settings &settings, const LArDriftVolumeList &driftVolumeList)
	Create pandora LArTPCs to represent the different drift volumes in use. More...

static void	CreatePandoraDetectorGaps (const Settings &settings, const LArDriftVolumeList &driftVolumeList, const LArDetectorGapList &listOfGaps)
	Create pandora line gaps to cover dead regions between TPCs in a global drift volume approach. More...

static void	CreatePandoraReadoutGaps (const Settings &settings, const LArDriftVolumeMap &driftVolumeMap)
	Create pandora line gaps to cover any (continuous regions of) bad channels. More...

static void	CreatePandoraMCParticles (const Settings &settings, const MCTruthToMCParticles &truthToParticles, const MCParticlesToMCTruth &particlesToTruth, const RawMCParticleVector &generatorMCParticleVector)
	Create the Pandora MC particles from the MC particles. More...

static void	FindPrimaryParticles (const RawMCParticleVector &mcParticleVector, std::map< const simb::MCParticle, bool > &primaryMCParticleMap)
	Find all primary MCParticles in a given vector of MCParticles. More...

static bool	IsPrimaryMCParticle (const art::Ptr< simb::MCParticle > &mcParticle, std::map< const simb::MCParticle, bool > &primaryMCParticleMap)
	Check whether an MCParticle can be found in a given map. More...

static void	CreatePandoraMCLinks2D (const Settings &settings, const HitMap &hitMap, const HitsToTrackIDEs &hitToParticleMap)
	Create links between the 2D hits and Pandora MC particles. More...

Private Types
typedef std::map< std::string, lar_content::MCProcess >	MCProcessMap

Static Private Member Functions
static void	GetTrueStartAndEndPoints (const Settings &settings, const art::Ptr< simb::MCParticle > &particle, int &startT, int &endT)
	Loop over MC trajectory points and identify start and end points within the detector. More...

static void	GetTrueStartAndEndPoints (const unsigned int cstat, const unsigned int tpc, const art::Ptr< simb::MCParticle > &particle, int &startT, int &endT)
	Loop over MC trajectory points and identify start and end points within a given cryostat and TPC. More...

static float	GetTrueX0 (const art::Event &evt, const art::Ptr< simb::MCParticle > &particle, const int nT)
	Use detector and time services to get a true X offset for a given trajectory point. More...

static double	GetMips (const detinfo::DetectorPropertiesData &detProp, const Settings &settings, const double hit_Charge, const geo::View_t hit_View)
	Convert charge in ADCs to approximate MIPs. More...

static void	FillMCProcessMap (MCProcessMap &processMap)
	Populate a map from MC process string to enumeration. More...

Detailed Description

LArPandoraInput class.

Definition at line 27 of file LArPandoraInput.h.

Member Typedef Documentation

typedef std::map<std::string, lar_content::MCProcess> lar_pandora::LArPandoraInput::MCProcessMap

private

Definition at line 141 of file LArPandoraInput.h.

Member Function Documentation

void lar_pandora::LArPandoraInput::CreatePandoraDetectorGaps	(	const Settings &	settings,
		const LArDriftVolumeList &	driftVolumeList,
		const LArDetectorGapList &	listOfGaps
	)

static

Create pandora line gaps to cover dead regions between TPCs in a global drift volume approach.

Parameters

settings	the settings
driftVolumeList	the drift volume list
listOfGaps	the list of gaps

Definition at line 245 of file LArPandoraInput.cxx.

   {
     //ATTN - Unlike SP, DP detector gaps are not in the drift direction
     art::ServiceHandle<geo::Geometry const> theGeometry;
     LArPandoraDetectorType* detType(detector_functions::GetDetectorType());
 
     mf::LogDebug("LArPandora") << " *** LArPandoraInput::CreatePandoraDetectorGaps(...) *** "
                                << std::endl;
 
     if (!settings.m_pPrimaryPandora)
       throw cet::exception("LArPandora")
         << "CreatePandoraDetectorGaps - primary Pandora instance does not exist ";
 
     const pandora::Pandora* pPandora(settings.m_pPrimaryPandora);
 
     for (const LArDetectorGap& gap : listOfGaps) {
       PandoraApi::Geometry::LineGap::Parameters parameters;
 
       try {
         parameters = detType->CreateLineGapParametersFromDetectorGaps(gap);
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora")
           << "CreatePandoraDetectorGaps - invalid line gap parameter provided, all assigned values "
              "must be finite, line gap omitted "
           << std::endl;
         continue;
       }
       try {
         PANDORA_THROW_RESULT_IF(pandora::STATUS_CODE_SUCCESS,
                                 !=,
                                 PandoraApi::Geometry::LineGap::Create(*pPandora, parameters));
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora") << "CreatePandoraDetectorGaps - unable to create line gap, "
                                         "insufficient or invalid information supplied "
                                      << std::endl;
         continue;
       }
     }
   }

void lar_pandora::LArPandoraInput::CreatePandoraHits2D	(	const art::Event &	evt,
		const Settings &	settings,
		const LArDriftVolumeMap &	driftVolumeMap,
		const HitVector &	hitVector,
		IdToHitMap &	idToHitMap
	)

static

Create the Pandora 2D hits from the ART hits.

Parameters

evt	art event being processed
settings	the settings
driftVolumeMap	the mapping from volume id to drift volume
hits	the input list of ART hits for this event
idToHitMap	to receive the mapping from Pandora hit ID to ART hit

Definition at line 45 of file LArPandoraInput.cxx.

   {
     mf::LogDebug("LArPandora") << " *** LArPandoraInput::CreatePandoraHits2D(...) *** "
                                << std::endl;
 
     if (!settings.m_pPrimaryPandora)
       throw cet::exception("LArPandora")
         << "CreatePandoraHits2D - primary Pandora instance does not exist ";
 
     const pandora::Pandora* pPandora(settings.m_pPrimaryPandora);
 
     art::ServiceHandle<geo::Geometry const> theGeometry;
     auto const detProp = art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(e);
     LArPandoraDetectorType* detType(detector_functions::GetDetectorType());
 
     // Loop over ART hits
     int hitCounter(settings.m_hitCounterOffset);
 
     lar_content::LArCaloHitFactory caloHitFactory;
 
     for (HitVector::const_iterator iter = hitVector.begin(), iterEnd = hitVector.end();
          iter != iterEnd;
          ++iter) {
       const art::Ptr<recob::Hit> hit = *iter;
       const geo::WireID hit_WireID(hit->WireID());
 
       // Get basic hit properties (view, time, charge)
       const geo::View_t hit_View(hit->View());
       const double hit_Charge(hit->Integral());
       const double hit_Time(hit->PeakTime());
       const double hit_TimeStart(hit->PeakTimeMinusRMS());
       const double hit_TimeEnd(hit->PeakTimePlusRMS());
 
       // Get hit X coordinate and, if using a single global drift volume, remove any out-of-time hits here
       const double xpos_cm(
         detProp.ConvertTicksToX(hit_Time, hit_WireID.Plane, hit_WireID.TPC, hit_WireID.Cryostat));
       const double dxpos_cm(
         std::fabs(detProp.ConvertTicksToX(
                     hit_TimeEnd, hit_WireID.Plane, hit_WireID.TPC, hit_WireID.Cryostat) -
                   detProp.ConvertTicksToX(
                     hit_TimeStart, hit_WireID.Plane, hit_WireID.TPC, hit_WireID.Cryostat)));
 
       // Get hit Y and Z coordinates, based on central position of wire
       double xyz[3];
       theGeometry->Cryostat(hit_WireID.Cryostat)
         .TPC(hit_WireID.TPC)
         .Plane(hit_WireID.Plane)
         .Wire(hit_WireID.Wire)
         .GetCenter(xyz);
       const double y0_cm(xyz[1]);
       const double z0_cm(xyz[2]);
 
       // Get other hit properties here
       const double wire_pitch_cm(theGeometry->WirePitch(hit_View)); // cm
       const double mips(LArPandoraInput::GetMips(detProp, settings, hit_Charge, hit_View));
 
       // Create Pandora CaloHit
       lar_content::LArCaloHitParameters caloHitParameters;
 
       try {
         caloHitParameters.m_expectedDirection = pandora::CartesianVector(0., 0., 1.);
         caloHitParameters.m_cellNormalVector = pandora::CartesianVector(0., 0., 1.);
         caloHitParameters.m_cellSize0 = settings.m_dx_cm;
         caloHitParameters.m_cellSize1 = (settings.m_useHitWidths ? dxpos_cm : settings.m_dx_cm);
         caloHitParameters.m_cellThickness = wire_pitch_cm;
         caloHitParameters.m_cellGeometry = pandora::RECTANGULAR;
         caloHitParameters.m_time = 0.;
         caloHitParameters.m_nCellRadiationLengths = settings.m_dx_cm / settings.m_rad_cm;
         caloHitParameters.m_nCellInteractionLengths = settings.m_dx_cm / settings.m_int_cm;
         caloHitParameters.m_isDigital = false;
         caloHitParameters.m_hitRegion = pandora::SINGLE_REGION;
         caloHitParameters.m_layer = 0;
         caloHitParameters.m_isInOuterSamplingLayer = false;
         caloHitParameters.m_inputEnergy = hit_Charge;
         caloHitParameters.m_mipEquivalentEnergy = mips;
         caloHitParameters.m_electromagneticEnergy = mips * settings.m_mips_to_gev;
         caloHitParameters.m_hadronicEnergy = mips * settings.m_mips_to_gev;
         caloHitParameters.m_pParentAddress = (void*)((intptr_t)(++hitCounter));
         caloHitParameters.m_larTPCVolumeId =
           LArPandoraGeometry::GetVolumeID(driftVolumeMap, hit_WireID.Cryostat, hit_WireID.TPC);
         caloHitParameters.m_daughterVolumeId = LArPandoraGeometry::GetDaughterVolumeID(
           driftVolumeMap, hit_WireID.Cryostat, hit_WireID.TPC);
 
         if (hit_View == detType->TargetViewW(hit_WireID.TPC, hit_WireID.Cryostat)) {
           caloHitParameters.m_hitType = pandora::TPC_VIEW_W;
           const double wpos_cm(
             pPandora->GetPlugins()->GetLArTransformationPlugin()->YZtoW(y0_cm, z0_cm));
           caloHitParameters.m_positionVector = pandora::CartesianVector(xpos_cm, 0., wpos_cm);
         }
         else if (hit_View == detType->TargetViewU(hit_WireID.TPC, hit_WireID.Cryostat)) {
           caloHitParameters.m_hitType = pandora::TPC_VIEW_U;
           const double upos_cm(
             pPandora->GetPlugins()->GetLArTransformationPlugin()->YZtoU(y0_cm, z0_cm));
           caloHitParameters.m_positionVector = pandora::CartesianVector(xpos_cm, 0., upos_cm);
         }
         else if (hit_View == detType->TargetViewV(hit_WireID.TPC, hit_WireID.Cryostat)) {
           caloHitParameters.m_hitType = pandora::TPC_VIEW_V;
           const double vpos_cm(
             pPandora->GetPlugins()->GetLArTransformationPlugin()->YZtoV(y0_cm, z0_cm));
           caloHitParameters.m_positionVector = pandora::CartesianVector(xpos_cm, 0., vpos_cm);
         }
         else {
           throw cet::exception("LArPandora")
             << "CreatePandoraHits2D - this wire view not recognised (View=" << hit_View << ") ";
         }
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora")
           << "CreatePandoraHits2D - invalid calo hit parameter provided, all assigned values must "
              "be finite, calo hit omitted "
           << std::endl;
         continue;
       }
 
       // Store the hit address
       if (hitCounter >= settings.m_uidOffset)
         throw cet::exception("LArPandora")
           << "CreatePandoraHits2D - detected an excessive number of hits (" << hitCounter << ") ";
 
       idToHitMap[hitCounter] = hit;
 
       // Create the Pandora hit
       try {
         PANDORA_THROW_RESULT_IF(
           pandora::STATUS_CODE_SUCCESS,
           !=,
           PandoraApi::CaloHit::Create(*pPandora, caloHitParameters, caloHitFactory));
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora") << "CreatePandoraHits2D - unable to create calo hit, "
                                         "insufficient or invalid information supplied "
                                      << std::endl;
         continue;
       }
     }
   }

void lar_pandora::LArPandoraInput::CreatePandoraLArTPCs	(	const Settings &	settings,
		const LArDriftVolumeList &	driftVolumeList
	)

static

Create pandora LArTPCs to represent the different drift volumes in use.

Parameters

settings	the settings
driftVolumeList	the drift volume list

Definition at line 189 of file LArPandoraInput.cxx.

   {
     mf::LogDebug("LArPandora") << " *** LArPandoraInput::CreatePandoraLArTPCs(...) *** "
                                << std::endl;
 
     if (!settings.m_pPrimaryPandora)
       throw cet::exception("LArPandora")
         << "CreatePandoraDetectorGaps - primary Pandora instance does not exist ";
 
     const pandora::Pandora* pPandora(settings.m_pPrimaryPandora);
 
     for (const LArDriftVolume& driftVolume : driftVolumeList) {
       PandoraApi::Geometry::LArTPC::Parameters parameters;
 
       try {
         parameters.m_larTPCVolumeId = driftVolume.GetVolumeID();
         parameters.m_centerX = driftVolume.GetCenterX();
         parameters.m_centerY = driftVolume.GetCenterY();
         parameters.m_centerZ = driftVolume.GetCenterZ();
         parameters.m_widthX = driftVolume.GetWidthX();
         parameters.m_widthY = driftVolume.GetWidthY();
         parameters.m_widthZ = driftVolume.GetWidthZ();
         parameters.m_wirePitchU = driftVolume.GetWirePitchU();
         parameters.m_wirePitchV = driftVolume.GetWirePitchV();
         parameters.m_wirePitchW = driftVolume.GetWirePitchW();
         parameters.m_wireAngleU = driftVolume.GetWireAngleU();
         parameters.m_wireAngleV = driftVolume.GetWireAngleV();
         parameters.m_wireAngleW = driftVolume.GetWireAngleW();
         parameters.m_sigmaUVW = driftVolume.GetSigmaUVZ();
         parameters.m_isDriftInPositiveX = driftVolume.IsPositiveDrift();
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora") << "CreatePandoraLArTPCs - invalid tpc parameter provided, "
                                         "all assigned values must be finite, tpc omitted "
                                      << std::endl;
         continue;
       }
 
       try {
         PANDORA_THROW_RESULT_IF(pandora::STATUS_CODE_SUCCESS,
                                 !=,
                                 PandoraApi::Geometry::LArTPC::Create(*pPandora, parameters));
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora") << "CreatePandoraLArTPCs - unable to create tpc, insufficient "
                                         "or invalid information supplied "
                                      << std::endl;
         continue;
       }
     }
   }

void lar_pandora::LArPandoraInput::CreatePandoraMCLinks2D	(	const Settings &	settings,
		const HitMap &	hitMap,
		const HitsToTrackIDEs &	hitToParticleMap
	)

static

Create links between the 2D hits and Pandora MC particles.

Parameters

settings	the settings
hitMap	mapping from Pandora hit addresses to ART hits
hitToParticleMap	mapping from each ART hit to its underlying G4 track ID

Definition at line 686 of file LArPandoraInput.cxx.

   {
     mf::LogDebug("LArPandora") << " *** LArPandoraInput::CreatePandoraMCLinks(...) *** "
                                << std::endl;
 
     if (!settings.m_pPrimaryPandora)
       throw cet::exception("LArPandora")
         << "CreatePandoraMCLinks2D - primary Pandora instance does not exist ";
 
     const pandora::Pandora* pPandora(settings.m_pPrimaryPandora);
 
     for (IdToHitMap::const_iterator iterI = idToHitMap.begin(), iterEndI = idToHitMap.end();
          iterI != iterEndI;
          ++iterI) {
       const int hitID(iterI->first);
       const art::Ptr<recob::Hit> hit(iterI->second);
       //  const geo::WireID hit_WireID(hit->WireID());
 
       // Get list of associated MC particles
       HitsToTrackIDEs::const_iterator iterJ = hitToParticleMap.find(hit);
 
       if (hitToParticleMap.end() == iterJ) continue;
 
       const TrackIDEVector& trackCollection = iterJ->second;
 
       if (trackCollection.size() == 0)
         throw cet::exception("LArPandora")
           << "CreatePandoraMCLinks2D - found a hit without any associated MC truth information ";
 
       // Create links between hits and MC particles
       for (unsigned int k = 0; k < trackCollection.size(); ++k) {
         const sim::TrackIDE trackIDE(trackCollection.at(k));
         const int trackID(std::abs(trackIDE.trackID)); // TODO: Find out why std::abs is needed
         const float energyFrac(trackIDE.energyFrac);
 
         try {
           PANDORA_THROW_RESULT_IF(
             pandora::STATUS_CODE_SUCCESS,
             !=,
             PandoraApi::SetCaloHitToMCParticleRelationship(
               *pPandora, (void*)((intptr_t)hitID), (void*)((intptr_t)trackID), energyFrac));
         }
         catch (const pandora::StatusCodeException&) {
           mf::LogWarning("LArPandora") << "CreatePandoraMCLinks2D - unable to create calo hit to "
                                           "mc particle relationship, invalid information supplied "
                                        << std::endl;
           continue;
         }
       }
     }
   }

void lar_pandora::LArPandoraInput::CreatePandoraMCParticles	(	const Settings &	settings,
		const MCTruthToMCParticles &	truthToParticles,
		const MCParticlesToMCTruth &	particlesToTruth,
		const RawMCParticleVector &	generatorMCParticleVector
	)

static

Create the Pandora MC particles from the MC particles.

Parameters

settings	the settings
truthToParticles	mapping from MC truth to MC particles
particlesToTruth	mapping from MC particles to MC truth

Definition at line 397 of file LArPandoraInput.cxx.

   {
     mf::LogDebug("LArPandora") << " *** LArPandoraInput::CreatePandoraMCParticles(...) *** "
                                << std::endl;
     art::ServiceHandle<cheat::ParticleInventoryService const> particleInventoryService;
 
     if (!settings.m_pPrimaryPandora)
       throw cet::exception("LArPandora")
         << "CreatePandoraMCParticles - primary Pandora instance does not exist ";
 
     const pandora::Pandora* pPandora(settings.m_pPrimaryPandora);
 
     // Make indexed list of MC particles
     MCParticleMap particleMap;
 
     for (MCParticlesToMCTruth::const_iterator iter = particleToTruthMap.begin(),
                                               iterEnd = particleToTruthMap.end();
          iter != iterEnd;
          ++iter) {
       const art::Ptr<simb::MCParticle> particle = iter->first;
       particleMap[particle->TrackId()] = particle;
     }
 
     // Loop over MC truth objects
     int neutrinoCounter(0);
 
     lar_content::LArMCParticleFactory mcParticleFactory;
 
     for (MCTruthToMCParticles::const_iterator iter1 = truthToParticleMap.begin(),
                                               iterEnd1 = truthToParticleMap.end();
          iter1 != iterEnd1;
          ++iter1) {
       const art::Ptr<simb::MCTruth> truth = iter1->first;
 
       if (truth->NeutrinoSet()) {
         const simb::MCNeutrino neutrino(truth->GetNeutrino());
         ++neutrinoCounter;
 
         if (neutrinoCounter >= settings.m_uidOffset)
           throw cet::exception("LArPandora")
             << "CreatePandoraMCParticles - detected an excessive number of mc neutrinos ("
             << neutrinoCounter << ")";
 
         const int neutrinoID(neutrinoCounter + 9 * settings.m_uidOffset);
 
         // Create Pandora 3D MC Particle
         lar_content::LArMCParticleParameters mcParticleParameters;
 
         try {
           mcParticleParameters.m_nuanceCode = neutrino.InteractionType();
           mcParticleParameters.m_process = lar_content::MC_PROC_INCIDENT_NU;
           mcParticleParameters.m_energy = neutrino.Nu().E();
           mcParticleParameters.m_momentum =
             pandora::CartesianVector(neutrino.Nu().Px(), neutrino.Nu().Py(), neutrino.Nu().Pz());
           mcParticleParameters.m_vertex =
             pandora::CartesianVector(neutrino.Nu().Vx(), neutrino.Nu().Vy(), neutrino.Nu().Vz());
           mcParticleParameters.m_endpoint =
             pandora::CartesianVector(neutrino.Nu().Vx(), neutrino.Nu().Vy(), neutrino.Nu().Vz());
           mcParticleParameters.m_particleId = neutrino.Nu().PdgCode();
           mcParticleParameters.m_mcParticleType = pandora::MC_3D;
           mcParticleParameters.m_pParentAddress = (void*)((intptr_t)neutrinoID);
         }
         catch (const pandora::StatusCodeException&) {
           mf::LogWarning("LArPandora")
             << "CreatePandoraMCParticles - invalid mc neutrino parameter provided, all assigned "
                "values must be finite, mc neutrino omitted "
             << std::endl;
           continue;
         }
 
         try {
           PANDORA_THROW_RESULT_IF(
             pandora::STATUS_CODE_SUCCESS,
             !=,
             PandoraApi::MCParticle::Create(*pPandora, mcParticleParameters, mcParticleFactory));
         }
         catch (const pandora::StatusCodeException&) {
           mf::LogWarning("LArPandora") << "CreatePandoraMCParticles - unable to create mc "
                                           "neutrino, insufficient or invalid information supplied "
                                        << std::endl;
           continue;
         }
 
         // Loop over associated particles
         const MCParticleVector& particleVector = iter1->second;
 
         for (MCParticleVector::const_iterator iter2 = particleVector.begin(),
                                               iterEnd2 = particleVector.end();
              iter2 != iterEnd2;
              ++iter2) {
           const art::Ptr<simb::MCParticle> particle = *iter2;
           const int trackID(particle->TrackId());
 
           // Mother/Daughter Links
           if (particle->Mother() == 0) {
             try {
               PANDORA_THROW_RESULT_IF(
                 pandora::STATUS_CODE_SUCCESS,
                 !=,
                 PandoraApi::SetMCParentDaughterRelationship(
                   *pPandora, (void*)((intptr_t)neutrinoID), (void*)((intptr_t)trackID)));
             }
             catch (const pandora::StatusCodeException&) {
               mf::LogWarning("LArPandora") << "CreatePandoraMCParticles - unable to create mc "
                                               "particle relationship, invalid information supplied "
                                            << std::endl;
               continue;
             }
           }
         }
       }
     }
 
     mf::LogDebug("LArPandora") << "   Number of Pandora neutrinos: " << neutrinoCounter
                                << std::endl;
 
     // Loop over G4 particles
     int particleCounter(0);
 
     // Find Primary Generator Particles
     std::map<const simb::MCParticle, bool> primaryGeneratorMCParticleMap;
     LArPandoraInput::FindPrimaryParticles(generatorMCParticleVector, primaryGeneratorMCParticleMap);
 
     for (MCParticleMap::const_iterator iterI = particleMap.begin(), iterEndI = particleMap.end();
          iterI != iterEndI;
          ++iterI) {
       const art::Ptr<simb::MCParticle> particle = iterI->second;
 
       if (particle->TrackId() != iterI->first)
         throw cet::exception("LArPandora") << "CreatePandoraMCParticles - mc truth information "
                                               "appears to be scrambled in this event";
 
       if (particle->TrackId() >= settings.m_uidOffset)
         throw cet::exception("LArPandora")
           << "CreatePandoraMCParticles - detected an excessive number of MC particles ("
           << particle->TrackId() << ")";
 
       ++particleCounter;
 
       // Find start and end trajectory points
       int firstT(-1), lastT(-1);
       LArPandoraInput::GetTrueStartAndEndPoints(settings, particle, firstT, lastT);
 
       if (firstT < 0 && lastT < 0) {
         firstT = 0;
         lastT = 0;
       }
 
       // Lookup position and kinematics at start and end points
       const float vtxX(particle->Vx(firstT));
       const float vtxY(particle->Vy(firstT));
       const float vtxZ(particle->Vz(firstT));
 
       const float endX(particle->Vx(lastT));
       const float endY(particle->Vy(lastT));
       const float endZ(particle->Vz(lastT));
 
       const float pX(particle->Px(firstT));
       const float pY(particle->Py(firstT));
       const float pZ(particle->Pz(firstT));
       const float E(particle->E(firstT));
 
       // Find the source of the mc particle
       int nuanceCode(0);
       const int trackID(particle->TrackId());
       const simb::Origin_t origin(particleInventoryService->TrackIdToMCTruth(trackID).Origin());
 
       if (LArPandoraInput::IsPrimaryMCParticle(particle, primaryGeneratorMCParticleMap)) {
         nuanceCode = 2001;
       }
       else if (simb::kCosmicRay == origin) {
         nuanceCode = 3000;
       }
       else if (simb::kSingleParticle == origin) {
         nuanceCode = 2000;
       }
 
       // Create 3D Pandora MC Particle
       lar_content::LArMCParticleParameters mcParticleParameters;
 
       try {
         MCProcessMap processMap;
         FillMCProcessMap(processMap);
         mcParticleParameters.m_nuanceCode = nuanceCode;
         if (processMap.find(particle->Process()) != processMap.end()) {
           mcParticleParameters.m_process = processMap[particle->Process()];
         }
         else {
           mcParticleParameters.m_process = lar_content::MC_PROC_UNKNOWN;
           mf::LogWarning("LArPandora")
             << "CreatePandoraMCParticles - found an unknown process" << std::endl;
         }
         mcParticleParameters.m_energy = E;
         mcParticleParameters.m_particleId = particle->PdgCode();
         mcParticleParameters.m_momentum = pandora::CartesianVector(pX, pY, pZ);
         mcParticleParameters.m_vertex = pandora::CartesianVector(vtxX, vtxY, vtxZ);
         mcParticleParameters.m_endpoint = pandora::CartesianVector(endX, endY, endZ);
         mcParticleParameters.m_mcParticleType = pandora::MC_3D;
         mcParticleParameters.m_pParentAddress = (void*)((intptr_t)particle->TrackId());
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora")
           << "CreatePandoraMCParticles - invalid mc particle parameter provided, all assigned "
              "values must be finite, mc particle omitted "
           << std::endl;
         continue;
       }
 
       try {
         PANDORA_THROW_RESULT_IF(
           pandora::STATUS_CODE_SUCCESS,
           !=,
           PandoraApi::MCParticle::Create(*pPandora, mcParticleParameters, mcParticleFactory));
       }
       catch (const pandora::StatusCodeException&) {
         mf::LogWarning("LArPandora") << "CreatePandoraMCParticles - unable to create mc particle, "
                                         "insufficient or invalid information supplied "
                                      << std::endl;
         continue;
       }
 
       // Create Mother/Daughter Links between 3D MC Particles
       const int id_mother(particle->Mother());
       MCParticleMap::const_iterator iterJ = particleMap.find(id_mother);
 
       if (iterJ != particleMap.end()) {
         try {
           PANDORA_THROW_RESULT_IF(
             pandora::STATUS_CODE_SUCCESS,
             !=,
             PandoraApi::SetMCParentDaughterRelationship(
               *pPandora, (void*)((intptr_t)id_mother), (void*)((intptr_t)particle->TrackId())));
         }
         catch (const pandora::StatusCodeException&) {
           mf::LogWarning("LArPandora") << "CreatePandoraMCParticles - Unable to create mc particle "
                                           "relationship, invalid information supplied "
                                        << std::endl;
           continue;
         }
       }
     }
 
     mf::LogDebug("LArPandora") << "Number of mc particles: " << particleCounter << std::endl;
   }

void lar_pandora::LArPandoraInput::CreatePandoraReadoutGaps	(	const Settings &	settings,
		const LArDriftVolumeMap &	driftVolumeMap
	)

static

Create pandora line gaps to cover any (continuous regions of) bad channels.

Parameters

settings	the settings
driftVolumeMap	the mapping from volume id to drift volume

Definition at line 292 of file LArPandoraInput.cxx.

   {
     mf::LogDebug("LArPandora") << " *** LArPandoraInput::CreatePandoraReadoutGaps(...) *** "
                                << std::endl;
 
     if (!settings.m_pPrimaryPandora)
       throw cet::exception("LArPandora")
         << "CreatePandoraReadoutGaps - primary Pandora instance does not exist ";
 
     const pandora::Pandora* pPandora(settings.m_pPrimaryPandora);
 
     art::ServiceHandle<geo::Geometry const> theGeometry;
     const lariov::ChannelStatusProvider& channelStatus(
       art::ServiceHandle<lariov::ChannelStatusService const>()->GetProvider());
 
     LArPandoraDetectorType* detType(detector_functions::GetDetectorType());
 
     for (unsigned int icstat = 0; icstat < theGeometry->Ncryostats(); ++icstat) {
       for (unsigned int itpc = 0; itpc < theGeometry->NTPC(icstat); ++itpc) {
         const geo::TPCGeo& TPC(theGeometry->TPC(itpc));
 
         for (unsigned int iplane = 0; iplane < TPC.Nplanes(); ++iplane) {
           const geo::PlaneGeo& plane(TPC.Plane(iplane));
           const float halfWirePitch(0.5f * theGeometry->WirePitch(plane.View()));
           const unsigned int nWires(theGeometry->Nwires(plane.ID()));
 
           int firstBadWire(-1), lastBadWire(-1);
 
           for (unsigned int iwire = 0; iwire < nWires; ++iwire) {
             const raw::ChannelID_t channel(
               theGeometry->PlaneWireToChannel(iplane, iwire, itpc, icstat));
             const bool isBadChannel(channelStatus.IsBad(channel));
             const bool isLastWire(nWires == (iwire + 1));
 
             if (isBadChannel && (firstBadWire < 0)) firstBadWire = iwire;
 
             if (isBadChannel || isLastWire) lastBadWire = iwire;
 
             if (isBadChannel && !isLastWire) continue;
 
             if ((firstBadWire < 0) || (lastBadWire < 0)) continue;
 
             double firstXYZ[3], lastXYZ[3];
             theGeometry->Cryostat(icstat)
               .TPC(itpc)
               .Plane(iplane)
               .Wire(firstBadWire)
               .GetCenter(firstXYZ);
             theGeometry->Cryostat(icstat)
               .TPC(itpc)
               .Plane(iplane)
               .Wire(lastBadWire)
               .GetCenter(lastXYZ);
 
             firstBadWire = -1;
             lastBadWire = -1;
 
             PandoraApi::Geometry::LineGap::Parameters parameters;
 
             try {
               float xFirst(-std::numeric_limits<float>::max());
               float xLast(std::numeric_limits<float>::max());
 
               const unsigned int volumeId(
                 LArPandoraGeometry::GetVolumeID(driftVolumeMap, icstat, itpc));
               LArDriftVolumeMap::const_iterator volumeIter(driftVolumeMap.find(volumeId));
 
               if (driftVolumeMap.end() != volumeIter) {
                 xFirst = volumeIter->second.GetCenterX() - 0.5f * volumeIter->second.GetWidthX();
                 xLast = volumeIter->second.GetCenterX() + 0.5f * volumeIter->second.GetWidthX();
               }
 
               const geo::View_t iview = plane.View();
               parameters = detType->CreateLineGapParametersFromReadoutGaps(
                 iview, itpc, icstat, firstXYZ, lastXYZ, halfWirePitch, xFirst, xLast, pPandora);
             }
             catch (const pandora::StatusCodeException&) {
               mf::LogWarning("LArPandora")
                 << "CreatePandoraReadoutGaps - invalid line gap parameter provided, all assigned "
                    "values must be finite, line gap omitted "
                 << std::endl;
               continue;
             }
 
             try {
               PANDORA_THROW_RESULT_IF(pandora::STATUS_CODE_SUCCESS,
                                       !=,
                                       PandoraApi::Geometry::LineGap::Create(*pPandora, parameters));
             }
             catch (const pandora::StatusCodeException&) {
               mf::LogWarning("LArPandora") << "CreatePandoraReadoutGaps - unable to create line "
                                               "gap, insufficient or invalid information supplied "
                                            << std::endl;
               continue;
             }
           }
         }
       }
     }
   }

void lar_pandora::LArPandoraInput::FillMCProcessMap ( MCProcessMap & processMap )

staticprivate

Populate a map from MC process string to enumeration.

Parameters

processMap the output map from MC process string to enumeration

Definition at line 852 of file LArPandoraInput.cxx.

   {
     // QGSP_BERT and EM standard physics list mappings
     processMap["unknown"] = lar_content::MC_PROC_UNKNOWN;
     processMap["primary"] = lar_content::MC_PROC_PRIMARY;
     processMap["compt"] = lar_content::MC_PROC_COMPT;
     processMap["phot"] = lar_content::MC_PROC_PHOT;
     processMap["annihil"] = lar_content::MC_PROC_ANNIHIL;
     processMap["eIoni"] = lar_content::MC_PROC_E_IONI;
     processMap["eBrem"] = lar_content::MC_PROC_E_BREM;
     processMap["conv"] = lar_content::MC_PROC_CONV;
     processMap["muIoni"] = lar_content::MC_PROC_MU_IONI;
     processMap["muMinusCaptureAtRest"] = lar_content::MC_PROC_MU_MINUS_CAPTURE_AT_REST;
     processMap["neutronInelastic"] = lar_content::MC_PROC_NEUTRON_INELASTIC;
     processMap["nCapture"] = lar_content::MC_PROC_N_CAPTURE;
     processMap["hadElastic"] = lar_content::MC_PROC_HAD_ELASTIC;
     processMap["Decay"] = lar_content::MC_PROC_DECAY;
     processMap["CoulombScat"] = lar_content::MC_PROC_COULOMB_SCAT;
     processMap["muBrems"] = lar_content::MC_PROC_MU_BREM;
     processMap["muPairProd"] = lar_content::MC_PROC_MU_PAIR_PROD;
     processMap["PhotonInelastic"] = lar_content::MC_PROC_PHOTON_INELASTIC;
     processMap["hIoni"] = lar_content::MC_PROC_HAD_IONI;
     processMap["protonInelastic"] = lar_content::MC_PROC_PROTON_INELASTIC;
     processMap["pi+Inelastic"] = lar_content::MC_PROC_PI_PLUS_INELASTIC;
     processMap["CHIPSNuclearCaptureAtRest"] = lar_content::MC_PROC_CHIPS_NUCLEAR_CAPTURE_AT_REST;
     processMap["pi-Inelastic"] = lar_content::MC_PROC_PI_MINUS_INELASTIC;
     processMap["Transportation"] = lar_content::MC_PROC_TRANSPORTATION;
     processMap["Rayl"] = lar_content::MC_PROC_RAYLEIGH;
     processMap["hBrems"] = lar_content::MC_PROC_HAD_BREM;
     processMap["hPairProd"] = lar_content::MC_PROC_HAD_PAIR_PROD;
     processMap["ionIoni"] = lar_content::MC_PROC_ION_IONI;
     processMap["nKiller"] = lar_content::MC_PROC_NEUTRON_KILLER;
     processMap["ionInelastic"] = lar_content::MC_PROC_ION_INELASTIC;
     processMap["He3Inelastic"] = lar_content::MC_PROC_HE3_INELASTIC;
     processMap["alphaInelastic"] = lar_content::MC_PROC_ALPHA_INELASTIC;
     processMap["anti_He3Inelastic"] = lar_content::MC_PROC_ANTI_HE3_INELASTIC;
     processMap["anti_alphaInelastic"] = lar_content::MC_PROC_ANTI_ALPHA_INELASTIC;
     processMap["hFritiofCaptureAtRest"] = lar_content::MC_PROC_HAD_FRITIOF_CAPTURE_AT_REST;
     processMap["anti_deuteronInelastic"] = lar_content::MC_PROC_ANTI_DEUTERON_INELASTIC;
     processMap["anti_neutronInelastic"] = lar_content::MC_PROC_ANTI_NEUTRON_INELASTIC;
     processMap["anti_protonInelastic"] = lar_content::MC_PROC_ANTI_PROTON_INELASTIC;
     processMap["anti_tritonInelastic"] = lar_content::MC_PROC_ANTI_TRITON_INELASTIC;
     processMap["dInelastic"] = lar_content::MC_PROC_DEUTERON_INELASTIC;
     processMap["electronNuclear"] = lar_content::MC_PROC_ELECTRON_NUCLEAR;
     processMap["photonNuclear"] = lar_content::MC_PROC_PHOTON_NUCLEAR;
     processMap["kaon+Inelastic"] = lar_content::MC_PROC_KAON_PLUS_INELASTIC;
     processMap["kaon-Inelastic"] = lar_content::MC_PROC_KAON_MINUS_INELASTIC;
     processMap["hBertiniCaptureAtRest"] = lar_content::MC_PROC_HAD_BERTINI_CAPTURE_AT_REST;
     processMap["lambdaInelastic"] = lar_content::MC_PROC_LAMBDA_INELASTIC;
     processMap["muonNuclear"] = lar_content::MC_PROC_MU_NUCLEAR;
     processMap["tInelastic"] = lar_content::MC_PROC_TRITON_INELASTIC;
     processMap["primaryBackground"] = lar_content::MC_PROC_PRIMARY_BACKGROUND;
   }

void lar_pandora::LArPandoraInput::FindPrimaryParticles	(	const RawMCParticleVector &	mcParticleVector,
		std::map< const simb::MCParticle, bool > &	primaryMCParticleMap
	)

static

Find all primary MCParticles in a given vector of MCParticles.

Parameters

mcParticleVector	vector of all MCParticles to consider
primaryMCParticleMap	map containing primary MCParticles and bool indicating whether particle has been accounted for

Definition at line 648 of file LArPandoraInput.cxx.

   {
     for (const simb::MCParticle& mcParticle : mcParticleVector) {
       if ("primary" == mcParticle.Process()) {
         primaryMCParticleMap.emplace(std::make_pair(mcParticle, false));
       }
     }
   }

double lar_pandora::LArPandoraInput::GetMips	(	const detinfo::DetectorPropertiesData &	detProp,
		const Settings &	settings,
		const double	hit_Charge,
		const geo::View_t	hit_View
	)

staticprivate

Convert charge in ADCs to approximate MIPs.

Parameters

settings	the settings
hit_Charge	the input charge
hit_View	the input view number

Definition at line 826 of file LArPandoraInput.cxx.

   {
     art::ServiceHandle<geo::Geometry const> theGeometry;
 
     // TODO: Unite this procedure with other calorimetry procedures under development
     const double dQdX(hit_Charge / (theGeometry->WirePitch(hit_View))); // ADC/cm
     const double dQdX_e(dQdX /
                         (detProp.ElectronsToADC() * settings.m_recombination_factor)); // e/cm
     const double dEdX(settings.m_useBirksCorrection ?
                         detProp.BirksCorrection(dQdX_e) :
                         dQdX_e * 1000. / util::kGeVToElectrons); // MeV/cm
     double mips(dEdX / settings.m_dEdX_mip);
 
     if (mips < 0.) mips = settings.m_mips_if_negative;
 
     if (mips > settings.m_mips_max) mips = settings.m_mips_max;
 
     return mips;
   }

void lar_pandora::LArPandoraInput::GetTrueStartAndEndPoints	(	const Settings &	settings,
		const art::Ptr< simb::MCParticle > &	particle,
		int &	startT,
		int &	endT
	)

staticprivate

Loop over MC trajectory points and identify start and end points within the detector.

Parameters

settings	the settings
particle	the true particle
startT	the first trajectory point in the detector
endT	the last trajectory point in the detector

Definition at line 743 of file LArPandoraInput.cxx.

   {
     art::ServiceHandle<geo::Geometry const> theGeometry;
     firstT = -1;
     lastT = -1;
 
     for (unsigned int icstat = 0; icstat < theGeometry->Ncryostats(); ++icstat) {
       for (unsigned int itpc = 0; itpc < theGeometry->NTPC(icstat); ++itpc) {
         int thisfirstT(-1), thislastT(-1);
         LArPandoraInput::GetTrueStartAndEndPoints(icstat, itpc, particle, thisfirstT, thislastT);
 
         if (thisfirstT < 0) continue;
 
         if (firstT < 0 || thisfirstT < firstT) firstT = thisfirstT;
 
         if (lastT < 0 || thislastT > lastT) lastT = thislastT;
       }
     }
   }

void lar_pandora::LArPandoraInput::GetTrueStartAndEndPoints	(	const unsigned int	cstat,
		const unsigned int	tpc,
		const art::Ptr< simb::MCParticle > &	particle,
		int &	startT,
		int &	endT
	)

staticprivate

Loop over MC trajectory points and identify start and end points within a given cryostat and TPC.

Parameters

cstat	the cryostat
tpc	the TPC
particle	the true particle
startT	the first trajectory point in the detector
endT	the last trajectory point in the detector

Definition at line 769 of file LArPandoraInput.cxx.

   {
     art::ServiceHandle<geo::Geometry const> theGeometry;
 
     bool foundStartPosition(false);
     const int numTrajectoryPoints(static_cast<int>(particle->NumberTrajectoryPoints()));
 
     for (int nt = 0; nt < numTrajectoryPoints; ++nt) {
       const double pos[3] = {particle->Vx(nt), particle->Vy(nt), particle->Vz(nt)};
       geo::TPCID tpcID = theGeometry->FindTPCAtPosition(pos);
 
       if (!tpcID.isValid) continue;
 
       if (!(cstat == tpcID.Cryostat && tpc == tpcID.TPC)) continue;
 
       endT = nt;
 
       if (!foundStartPosition) {
         startT = endT;
         foundStartPosition = true;
       }
     }
   }

float lar_pandora::LArPandoraInput::GetTrueX0	(	const art::Event &	evt,
		const art::Ptr< simb::MCParticle > &	particle,
		const int	nT
	)

staticprivate

Use detector and time services to get a true X offset for a given trajectory point.

Parameters

evt	event currently being processing by art
particle	the true particle
nT	the trajectory point

Definition at line 800 of file LArPandoraInput.cxx.

   {
     art::ServiceHandle<geo::Geometry const> theGeometry;
     auto const clock_data = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(e);
     auto const det_prop =
       art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(e, clock_data);
 
     unsigned int which_tpc(0);
     unsigned int which_cstat(0);
     double pos[3] = {particle->Vx(nt), particle->Vy(nt), particle->Vz(nt)};
     theGeometry->PositionToTPC(pos, which_tpc, which_cstat);
 
     const float vtxT(particle->T(nt));
     const float vtxTDC(clock_data.TPCG4Time2Tick(vtxT));
     const float vtxTDC0(trigger_offset(clock_data));
 
     const geo::TPCGeo& theTpc = theGeometry->Cryostat(which_cstat).TPC(which_tpc);
     const float driftDir((theTpc.DriftDirection() == geo::kNegX) ? +1.0 : -1.0);
     return (driftDir * (vtxTDC - vtxTDC0) * det_prop.GetXTicksCoefficient());
   }

bool lar_pandora::LArPandoraInput::IsPrimaryMCParticle	(	const art::Ptr< simb::MCParticle > &	mcParticle,
		std::map< const simb::MCParticle, bool > &	primaryMCParticleMap
	)

static

Check whether an MCParticle can be found in a given map.

Parameters

mcParticle	target MCParticle
primaryMCParticleMap	map containing primary MCParticles and bool indicating whether particle has been accounted for

Definition at line 662 of file LArPandoraInput.cxx.

   {
     for (auto& mcParticleIter : primaryMCParticleMap) {
       if (!mcParticleIter.second) {
         const simb::MCParticle primaryMCParticle(mcParticleIter.first);
 
         if (std::fabs(primaryMCParticle.Px() - mcParticle->Px()) <
               std::numeric_limits<double>::epsilon() &&
             std::fabs(primaryMCParticle.Py() - mcParticle->Py()) <
               std::numeric_limits<double>::epsilon() &&
             std::fabs(primaryMCParticle.Pz() - mcParticle->Pz()) <
               std::numeric_limits<double>::epsilon()) {
           mcParticleIter.second = true;
           return true;
         }
       }
     }
     return false;
   }

The documentation for this class was generated from the following files:

Classes

Static Public Member Functions

Private Types

Static Private Member Functions

Detailed Description

Member Typedef Documentation

Member Function Documentation