From 46f43ee7fc9fa878161651893003c2b0d96718e2 Mon Sep 17 00:00:00 2001 From: Yuri Oksuzian Date: Mon, 4 Nov 2019 21:54:33 -0600 Subject: [PATCH 1/2] Adding upper momentum cut in FilterStepPointMomentum to be used in cosmic resampling jobs --- .../src/StepPointMC1stHitDumper_module.cc | 26 ++- Filters/src/FilterStepPointMomentum_module.cc | 4 +- JobConfig/cosmic/dsstops_resampler.fcl | 19 +- JobConfig/mixing/DS-flash-TrkCal-cut.fcl | 1 + Mu2eG4/geom/geom_common_phase1.txt | 163 +----------------- 5 files changed, 40 insertions(+), 173 deletions(-) diff --git a/Analyses/src/StepPointMC1stHitDumper_module.cc b/Analyses/src/StepPointMC1stHitDumper_module.cc index 779b7d8a28..8d32e215ac 100644 --- a/Analyses/src/StepPointMC1stHitDumper_module.cc +++ b/Analyses/src/StepPointMC1stHitDumper_module.cc @@ -71,6 +71,11 @@ namespace mu2e { float x; float y; float z; + + float InitX; + float InitY; + float InitZ; + float time; float px; @@ -83,11 +88,17 @@ namespace mu2e { int pdgId; unsigned particleId; unsigned volumeCopyNumber; + int eventId; + int subrunId; VDHit() : x(std::numeric_limits::quiet_NaN()) , y(std::numeric_limits::quiet_NaN()) , z(std::numeric_limits::quiet_NaN()) + , InitX(std::numeric_limits::quiet_NaN()) + , InitY(std::numeric_limits::quiet_NaN()) + , InitZ(std::numeric_limits::quiet_NaN()) + , time(std::numeric_limits::quiet_NaN()) , px(std::numeric_limits::quiet_NaN()) @@ -100,14 +111,20 @@ namespace mu2e { , pdgId(0) , particleId(-1U) , volumeCopyNumber(-1U) + , eventId(0) + , subrunId(0) {} //---------------------------------------------------------------- - VDHit(const SimParticleTimeOffset& toff, const StepPointMC& hit) + VDHit(const SimParticleTimeOffset& toff, const art::Event& event, const StepPointMC& hit) : x(hit.position().x()) , y(hit.position().y()) , z(hit.position().z()) + , InitX(hit.simParticle()->startPosition().x()) + , InitY(hit.simParticle()->startPosition().y()) + , InitZ(hit.simParticle()->startPosition().z()) + , time(toff.timeWithOffsetsApplied(hit)) , px(hit.momentum().x()) @@ -122,6 +139,9 @@ namespace mu2e { , pdgId(hit.simParticle()->pdgId()) , particleId(hit.simParticle()->id().asUint()) , volumeCopyNumber(hit.volumeId()) + + , eventId(event.id().event()) + , subrunId(event.id().subRun()) {} }; // struct VDHit @@ -169,7 +189,7 @@ namespace mu2e { //================================================================ void StepPointMC1stHitDumper::beginJob() { art::ServiceHandle tfs; - static const char branchDesc[] = "x/F:y/F:z/F:time/F:px/F:py/F:pz/F:pmag/F:ek/F:charge/F:pdgId/I:particleId/i:volumeCopy/i"; + static const char branchDesc[] = "x/F:y/F:z/F:InitX/F:InitY/F:InitZ/F:time/F:px/F:py/F:pz/F:pmag/F:ek/F:charge/F:pdgId/I:particleId/i:volumeCopy/i:eventId/I:subrunId/I"; nt_ = tfs->make( "nt", "StepPointMC1stHitDumper ntuple"); nt_->Branch("hits", &hit_, branchDesc); if(writeProperTime_) { @@ -195,7 +215,7 @@ namespace mu2e { // std::cout << event.id() << ": " << stepPoints.at(i) << std::endl; // } if (stepPoints.size()) { - hit_ = VDHit(toff_, stepPoints.at(0)); + hit_ = VDHit(toff_, event, stepPoints.at(0)); // std::cout << event.id() << ", " << hit_.pdgId // << ", " << hit_.time << ", " << hit_.charge // << std::endl; diff --git a/Filters/src/FilterStepPointMomentum_module.cc b/Filters/src/FilterStepPointMomentum_module.cc index e6b95c3291..947455a3aa 100644 --- a/Filters/src/FilterStepPointMomentum_module.cc +++ b/Filters/src/FilterStepPointMomentum_module.cc @@ -23,6 +23,7 @@ namespace mu2e { typedef std::vector InputTags; InputTags inputs_; double cutMomentumMin_; + double cutMomentumMax_; // statistics counters unsigned numInputEvents_; @@ -37,6 +38,7 @@ namespace mu2e { FilterStepPointMomentum::FilterStepPointMomentum(const fhicl::ParameterSet& pset) : art::EDFilter{pset} , cutMomentumMin_(pset.get("cutMomentumMin")) + , cutMomentumMax_(pset.get("cutMomentumMax", std::numeric_limits::max())) , numInputEvents_(0) , numPassedEvents_(0) { @@ -53,7 +55,7 @@ namespace mu2e { for(const auto& cn : inputs_) { auto ih = event.getValidHandle(cn); for(const auto& hit : *ih) { - if(hit.momentum().mag() > cutMomentumMin_) { + if(hit.momentum().mag() > cutMomentumMin_ && hit.momentum().mag() < cutMomentumMax_) { passed = true; break; } diff --git a/JobConfig/cosmic/dsstops_resampler.fcl b/JobConfig/cosmic/dsstops_resampler.fcl index f62eb301b7..7dc06d8c42 100644 --- a/JobConfig/cosmic/dsstops_resampler.fcl +++ b/JobConfig/cosmic/dsstops_resampler.fcl @@ -109,13 +109,6 @@ physics.filters.dsFilter: { vetoDaughters : [] } -physics.filters.stepPointMomentumFilter: { - module_type: FilterStepPointMomentum - inputs : [ "dsResample:crvStage1" ] - cutMomentumMin: 45. // Drop low momentum tracks - cutMomentumMax: 500. // Drop straight tracks -} - physics.filters.trackerStepPointFilter: { module_type : TrackerStepPointFilter } @@ -128,6 +121,13 @@ physics.filters.detectorFilter: { vetoDaughters: [] } +physics.filters.stepPointMomentumFilter: { + module_type: FilterStepPointMomentum + inputs : [ "detectorFilter:tracker" ] + cutMomentumMin: 45. // Drop low momentum tracks + cutMomentumMax: 500. // Drop straight tracks +} + physics.filters.g4status: { module_type: FilterStatusG4 input: "g4run" @@ -135,7 +135,7 @@ physics.filters.g4status: { } -physics.trig: [genCounter, dsResample, g4run, g4consistent, stepPointMomentumFilter, dsFilter, trackerStepPointFilter, detectorFilter] +physics.trig: [genCounter, dsResample, g4run, g4consistent, dsFilter, trackerStepPointFilter, detectorFilter, stepPointMomentumFilter] physics.g4StatusFilter : [dsResample, g4run, "!g4status"] physics.outputs: [truncatedEvtsOutput, DSVacuumOut] physics.an: [genCountLogger] @@ -189,5 +189,4 @@ physics.producers.g4run.Mu2eG4CommonCut:{ physics.producers.g4run.ResourceLimits.maxSimParticleCollectionSize: 100000 // Test: //physics.filters.dsResample.fileNames: ["/pnfs/mu2e/persistent/users/oksuzian/workflow/cry_rs1_1019_g4_10_5/outstage/363672.fcllist_191019181515/00/00000/sim.oksuzian.minedep-filter-cat.cryconcat_rs1_1019.002701_00000000.art"] -//physics.filters.dsResample.fileNames: ["/pnfs/mu2e/scratch/users/oksuzian/workflow/cry_rs1_1019_g4_10_5/good/682690.fcllist_191024104002/00/00047/sim.oksuzian.minedep-filter.cry-minedep_filter-10-5.002701_00020094.art"] -//physics.filters.dsResample.fileNames: ["/pnfs/mu2e/scratch/users/oksuzian/workflow/cry_rs1_1019_g4_10_5/good/683306.fcllist_191024104850/00/00000/sim.oksuzian.minedep-filter.cry-minedep_filter-10-5.002701_00000385.art"] \ No newline at end of file +//physics.filters.dsResample.fileNames: ["/pnfs/mu2e/persistent/users/oksuzian/workflow/cry_rs1_1019_g4_10_5/outstage/683306.fcllist_191024104850/00/00000/sim.oksuzian.minedep-filter.cry-minedep_filter-10-5.002701_00000385.art"] \ No newline at end of file diff --git a/JobConfig/mixing/DS-flash-TrkCal-cut.fcl b/JobConfig/mixing/DS-flash-TrkCal-cut.fcl index a4f4dc3c53..007aca1a6e 100644 --- a/JobConfig/mixing/DS-flash-TrkCal-cut.fcl +++ b/JobConfig/mixing/DS-flash-TrkCal-cut.fcl @@ -51,3 +51,4 @@ services.TFileService.fileName : "nts.owner.DS-flash-TrkCal-cut.version.sequence services.message.destinations.log.categories.ArtReport.reportEvery : 1 services.message.destinations.log.categories.ArtReport.limit : 1 services.message.destinations.log.categories.ArtReport.timespan : 300 +services.GeometryService.simulatedDetector.tool_type : "Mu2e" diff --git a/Mu2eG4/geom/geom_common_phase1.txt b/Mu2eG4/geom/geom_common_phase1.txt index 95ef5eeb61..5ffdd4bc55 100644 --- a/Mu2eG4/geom/geom_common_phase1.txt +++ b/Mu2eG4/geom/geom_common_phase1.txt @@ -1,164 +1,9 @@ -// -// Top level geometry file first used in pass 2 of (stage 4 and after) of the CD3c simulation campagin. -// -// Warning: do not write 10000. as 10,000.; it will be read as two numbers (10., 0.). +// A geometry file for using all regular concrete shielding design. +// Most of the geometry is inherited from the existing shielding geometry, except as below. -string detector.name = "g4geom_v00"; - -bool hasHall = true; -bool hasTarget = true; -bool hasProtonAbsorber = true; -bool hasTSdA = true; -bool hasExternalShielding = true; -bool hasDiskCalorimeter = true; -bool hasBFieldManager = true; -bool hasBeamline = true; -bool hasVirtualDetector = true; // some components, e.g. ProtonAbsorber assume vd presence now; -bool hasCosmicRayShield = true; -bool hasSTM = true; -bool hasMBS = true; // note the two subcomponents, see mbs section below; - // no MBS implies no downstream hole in Cosmic Ray Passive Shield - // Magnetic field may be affected as well - -#include "Mu2eG4/geom/g4_visOptions.txt" - -//------------------------------------------- -// Mu2e geometry includes -//------------------------------------------- - -// X-offset of the PS(+x) and DS(-x) from the Mu2e origin. -// The origin of the detector coordinate system is on the DS axis at the specified z. -double mu2e.solenoidOffset = 3904.; // mm -double mu2e.detectorSystemZ0 = 10171.; // mm G4BL: (17730-7292=9801 mm) - -#include "Mu2eG4/geom/mu2eWorld.txt" -// mu2eHall.txt should be used with protonBeamDump_v02.txt, below -//#include "Mu2eG4/geom/mu2eHall.txt" -// whereas mu2eHall_v02.txt should be used with protonBeamDump_v03.txt, below -#include "Mu2eG4/geom/mu2eHall_v02.txt" - -// Solenoids -#include "Mu2eG4/geom/DetectorSolenoid_v03.txt" -#include "Mu2eG4/geom/DSShielding_v01.txt" -#include "Mu2eG4/geom/ProductionSolenoid_v02.txt" -#include "Mu2eG4/geom/psEnclosure_v04.txt" -#include "Mu2eG4/geom/PSShield_v06.txt" -#include "Mu2eG4/geom/PSExternalShielding_v01.txt" -#include "Mu2eG4/geom/TransportSolenoid_v03.txt" +// First, the base geometry to start from... +#include "Mu2eG4/geom/geom_common_current.txt" // External Shielding #include "Mu2eG4/geom/ExtShieldUpstream_v06.txt" #include "Mu2eG4/geom/ExtShieldDownstream_v06.txt" -#include "Mu2eG4/geom/Saddle_v02.txt" -#include "Mu2eG4/geom/Pipe_v03.txt" -#include "Mu2eG4/geom/ElectronicRack_v01.txt" - -//#include "Mu2eG4/geom/stoppingTarget_TDR.txt" // 17 foil tapered muon stopping target used for the TDR simulations -//#include "Mu2eG4/geom/stoppingTarget_CD3C_34foils.txt" // 34 foil muon stopping target to be used for the CDC3 simulations -#include "Mu2eG4/geom/stoppingTargetHoles_DOE_review_2017.txt" // 37 foil muon stopping target with holes - -#include "Mu2eG4/geom/TSdA_v01.txt" -#include "Mu2eG4/geom/muonBeamStop_v05.txt" - -//#include "Mu2eG4/geom/MSTM_v01.txt" // muon stopping target monitor (deprecated) -//#include "Mu2eG4/geom/STM_v01.txt" // (muon) stopping target monitor -//#include "Mu2eG4/geom/STM_v02.txt" // (muon) stopping target monitor -#include "Mu2eG4/geom/STM_v03.txt" // (muon) stopping target monitor - -// Proton Absorber -#include "Mu2eG4/geom/protonAbsorber_cylindrical_v02.txt" -#include "Mu2eG4/geom/degrader_v01.txt" - -#include "Mu2eG4/geom/ProductionTargetInPS.txt" -//#include "Mu2eG4/geom/protonBeamDump_v02.txt" -#include "Mu2eG4/geom/protonBeamDump_v03.txt" -#include "Mu2eG4/geom/extmon_fnal_v02.txt" - -#include "Mu2eG4/geom/tracker_v5.txt" - -// Crystal calorimeter -#include "Mu2eG4/geom/calorimeter_CsI.txt" - -//CRV counters -#include "Mu2eG4/geom/crv_counters_v07.txt" - - -//------------------------------------------- -// Magnetic field -//------------------------------------------- - -// Form of DS field: 0 is full field; -// 1 is full upstream, const downstream; -// 2 is const throughout -int detSolFieldForm = 0; -vector bfield.innerMaps = { - "BFieldMaps/Mau13/DSMap.header", - "BFieldMaps/Mau13/PSMap.header", - "BFieldMaps/Mau13/TSuMap.header", - "BFieldMaps/Mau13/TSdMap.header", - "BFieldMaps/Mau13/PStoDumpAreaMap.header", - "BFieldMaps/Mau13/ProtonDumpAreaMap.header", - "BFieldMaps/Mau13/DSExtension.header" -}; - -// Value of the uniform magnetic field with the DS volume (only for -// detSolFieldForm>0) -double toyDS.bz = 1.0; - -// Gradient of field in DS2 volume. Applied only in the case -// of detSolFieldForm=1 or detSolFieldForm=2. -double toyDS.gradient = 0.0; // Tesla/m - -// This is recommended field map. See geom_mecofield.txt to use the meco field. -string bfield.format = "G4BL"; - -// The other option is "meco" -string bfield.interpolationStyle = trilinear; - -int bfield.verbosityLevel = 0; -bool bfield.writeG4BLBinaries = false; - -vector bfield.outerMaps = { - "BFieldMaps/Mau9/ExtMonUCIInternal1AreaMap.header", - "BFieldMaps/Mau9/ExtMonUCIInternal2AreaMap.header", - "BFieldMaps/Mau9/ExtMonUCIAreaMap.header", - "BFieldMaps/Mau9/PSAreaMap.header" -}; - -// This scale factor is of limited use. -// It can make approximate sense to scale the PS field to get a rough -// answer; the answer will be wrong in detail. -// It never makes sense to scale the TS field. -// Not sure if it ever makes sense to scale the PS field. -double bfield.scaleFactor = 1.0; - -//--------------------------------------- -// Virtual detectors -//--------------------------------------- -double vd.halfLength = 0.01; //mm -int vd.verbosityLevel = 0; -bool vd.visible = true; -bool vd.solid = false; - -// // VD right in front of a hall wall -// double vd.ExtMonCommonPlane.z = -11999.99; - - -//--------------------------------------- -// Region visualization -//--------------------------------------- -#include "Mu2eG4/geom/visualization_regions.txt" - - -// -// -// End notes: -// -// 1) Sources of information: -// -// -// -// This tells emacs to view this file in c++ mode. -// Local Variables: -// mode:c++ -// End: From 1efa439f77783c76775ad2941e78a7c6af5ebceb Mon Sep 17 00:00:00 2001 From: Yuri Oksuzian Date: Tue, 5 Nov 2019 15:35:03 -0600 Subject: [PATCH 2/2] Changed to Config pattern from ParameterSet --- Filters/src/FilterStepPointMomentum_module.cc | 42 +++++++++++++++---- 1 file changed, 34 insertions(+), 8 deletions(-) diff --git a/Filters/src/FilterStepPointMomentum_module.cc b/Filters/src/FilterStepPointMomentum_module.cc index 947455a3aa..dec963e486 100644 --- a/Filters/src/FilterStepPointMomentum_module.cc +++ b/Filters/src/FilterStepPointMomentum_module.cc @@ -29,24 +29,50 @@ namespace mu2e { unsigned numInputEvents_; unsigned numPassedEvents_; public: - explicit FilterStepPointMomentum(const fhicl::ParameterSet& pset); + + struct Config { + using Name=fhicl::Name; + using Comment=fhicl::Comment; + + fhicl::Sequence inputs { + Name("inputs"), + Comment("Particles and StepPointMCs mentioned in thise collections will be preserved.") + }; + + fhicl::Atom cutMomentumMin { + Name("cutMomentumMin"), + Comment("The filter passes events if any of the step points satisties pmag>cutMomentumMin\n" + "By default cutMomentumMin=-inf\n"), + -std::numeric_limits::max() + }; + + fhicl::Atom cutMomentumMax { + Name("cutMomentumMax"), + Comment("The filter passes events if any of the step points satisties pmag::max() + }; + + }; + + using Parameters = art::EDFilter::Table; + explicit FilterStepPointMomentum(const Parameters& conf); virtual bool filter(art::Event& event) override; virtual void endJob() override; }; //================================================================ - FilterStepPointMomentum::FilterStepPointMomentum(const fhicl::ParameterSet& pset) - : art::EDFilter{pset} - , cutMomentumMin_(pset.get("cutMomentumMin")) - , cutMomentumMax_(pset.get("cutMomentumMax", std::numeric_limits::max())) + FilterStepPointMomentum::FilterStepPointMomentum(const Parameters& conf) + : art::EDFilter{conf} + , cutMomentumMin_(conf().cutMomentumMin()) + , cutMomentumMax_(conf().cutMomentumMax()) , numInputEvents_(0) , numPassedEvents_(0) { - typedef std::vector VS; - const VS in(pset.get("inputs")); - for(const auto& i : in) { + for(const auto& i : conf().inputs()) { inputs_.emplace_back(i); } + } //================================================================