Public Member Functions | Private Attributes | List of all members
AnalysisExample::dEdx Class Reference
Inheritance diagram for AnalysisExample::dEdx:
art::EDAnalyzer art::detail::Analyzer art::detail::LegacyModule art::Observer art::ModuleBase

Public Member Functions

 dEdx (fhicl::ParameterSet const &pset)
 
virtual ~dEdx ()
 
void beginJob ()
 
void beginRun (const art::Run &run)
 
void reconfigure (fhicl::ParameterSet const &pset)
 
void analyze (const art::Event &evt)
 
void endJob ()
 
void findDepositionVertex (std::vector< sim::SimChannel > SCHandle)
 
void checkDepositionVertex (std::vector< sim::SimChannel > SCHandle, double x, double y, double z)
 
- Public Member Functions inherited from art::EDAnalyzer
 EDAnalyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 EDAnalyzer (Table< Config > const &config)
 
std::string workerType () const
 
- Public Member Functions inherited from art::detail::Analyzer
virtual ~Analyzer () noexcept
 
 Analyzer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 Analyzer (Table< Config > const &config)
 
void doBeginJob (SharedResources const &resources)
 
void doEndJob ()
 
void doRespondToOpenInputFile (FileBlock const &fb)
 
void doRespondToCloseInputFile (FileBlock const &fb)
 
void doRespondToOpenOutputFiles (FileBlock const &fb)
 
void doRespondToCloseOutputFiles (FileBlock const &fb)
 
bool doBeginRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doEndRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)
 
- Public Member Functions inherited from art::Observer
 ~Observer () noexcept
 
 Observer (Observer const &)=delete
 
 Observer (Observer &&)=delete
 
Observeroperator= (Observer const &)=delete
 
Observeroperator= (Observer &&)=delete
 
void registerProducts (ProductDescriptions &, ModuleDescription const &)
 
void fillDescriptions (ModuleDescription const &)
 
fhicl::ParameterSetID selectorConfig () const
 
- Public Member Functions inherited from art::ModuleBase
virtual ~ModuleBase () noexcept
 
 ModuleBase ()
 
ModuleDescription const & moduleDescription () const
 
void setModuleDescription (ModuleDescription const &)
 
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables () const
 
void sortConsumables (std::string const &current_process_name)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Private Attributes

art::ServiceHandle< art::TFileService > tfs
 
art::ServiceHandle< geo::GeometryfGeom
 
std::string fSimulationProducerLabel
 
TH1D * fMCdEdxHist
 
std::vector< double > fMCdEdxVec
 
double fMaxMCdEdx = 0.
 
double dE
 
double trackLength
 
TH3D * fStartEdepHist
 
TH3D * fTrackEdepHist
 
TH3D * fShowerEdepHist
 
std::vector< unsigned int > fNullList
 
double xi = 0.
 
double yi = 0.
 
double zi = 1000.
 
double ziFalse = -500.
 
unsigned int nearDeps = 0
 
bool goodDepVertex = false
 
double fConvDist
 
TH1D * fConvDistHist
 
double Vxi
 
double Vyi
 
double Vzi
 

Additional Inherited Members

- Public Types inherited from art::EDAnalyzer
using WorkerType = WorkerT< EDAnalyzer >
 
using ModuleType = EDAnalyzer
 
- Protected Member Functions inherited from art::Observer
std::string const & processName () const
 
bool wantAllEvents () const noexcept
 
bool wantEvent (ScheduleID id, Event const &e) const
 
Handle< TriggerResultsgetTriggerResults (Event const &e) const
 
 Observer (fhicl::ParameterSet const &config)
 
 Observer (std::vector< std::string > const &select_paths, std::vector< std::string > const &reject_paths, fhicl::ParameterSet const &config)
 
- Protected Member Functions inherited from art::ModuleBase
ConsumesCollectorconsumesCollector ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > consumes (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 

Detailed Description

Definition at line 53 of file dEdx_module.cc.

Constructor & Destructor Documentation

AnalysisExample::dEdx::dEdx ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 112 of file dEdx_module.cc.

113  : EDAnalyzer(parameterSet)
114  {
115  this->reconfigure(parameterSet);
116  }
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:25
void reconfigure(fhicl::ParameterSet const &pset)
Definition: dEdx_module.cc:158
AnalysisExample::dEdx::~dEdx ( )
virtual

Definition at line 120 of file dEdx_module.cc.

121  {
122  }

Member Function Documentation

void AnalysisExample::dEdx::analyze ( const art::Event evt)

Definition at line 165 of file dEdx_module.cc.

166  {
167 
168  auto simChanHandle = event.getHandle<std::vector<sim::SimChannel>>(fSimulationProducerLabel);
169 
170  // initialize at the start of every event
171  dE = 0.;
172  trackLength = 0.;
173  goodDepVertex = false;
174  unsigned int nTries = 0;
175  ziFalse = -500.;
176 
177  // find conversion point
178  while(!goodDepVertex && nTries < 300)
179  {
180  this->findDepositionVertex( (*simChanHandle) );
181  this->checkDepositionVertex( (*simChanHandle), xi, yi, zi );
182  nTries++;
183  }
184 
185  std::cout << "Number of attemps at finding deposition vertex = " << nTries << "." << std::endl;
186 
187  // find energy deposited into first 2.5 cm of "track"
188  if(nTries < 300 && goodDepVertex)
189  {
190  fStartEdepHist->Fill(xi,yi,zi);
191  for(auto const& channel : (*simChanHandle))
192  {
193  if(fGeom->SignalType(channel.Channel()) == geo::kCollection)
194  {
195  auto const& timeSlices = channel.TDCIDEMap();
196  for(auto const& t : timeSlices)
197  {
198  auto const& eDeps = t.second;
199  for(auto const& eDep : eDeps)
200  {
201  fShowerEdepHist->Fill(eDep.x,eDep.y,eDep.z);
202  trackLength = std::sqrt( (eDep.x-xi)*(eDep.x-xi) + (eDep.y-yi)*(eDep.y-yi) + (eDep.z-zi)*(eDep.z-zi) );
203  if(trackLength < 2.5 /* && eDep.z >= zi*/)
204  {
205  dE += eDep.energy;
206  fTrackEdepHist->Fill(eDep.x,eDep.y,eDep.z);
207  }
208  } // every energy deposition
209  } // in every time slice
210  } // for the collection plane
211  } // and for every wire
212  dE /= 2.5;
213  fMCdEdxVec.push_back(dE);
214  if(dE > fMaxMCdEdx) fMaxMCdEdx = dE;
215  if(dE < 0.5) fNullList.push_back(event.id().event());
216  }
217 
218  auto particleHandle = event.getHandle< std::vector<simb::MCParticle> >(fSimulationProducerLabel);
219 
220  for(auto const& particle : (*particleHandle))
221  {
222  if(particle.Process() == "primary")
223  {
224  Vxi = particle.Vx();
225  Vyi = particle.Vy();
226  Vzi = particle.Vz();
227  }
228  }
229  fConvDist = std::sqrt( (Vxi-xi)*(Vxi-xi) + (Vyi-yi)*(Vyi-yi) + (Vzi-zi)*(Vzi-zi) );
230  fConvDistHist->Fill(fConvDist);
231  }
art::ServiceHandle< geo::Geometry > fGeom
Definition: dEdx_module.cc:70
std::string fSimulationProducerLabel
Definition: dEdx_module.cc:71
uint8_t channel
Definition: CRTFragment.hh:201
SigType_t SignalType(geo::PlaneID const &pid) const
Returns the type of signal on the channels of specified TPC plane.
void checkDepositionVertex(std::vector< sim::SimChannel > SCHandle, double x, double y, double z)
Definition: dEdx_module.cc:275
void findDepositionVertex(std::vector< sim::SimChannel > SCHandle)
Definition: dEdx_module.cc:249
std::vector< unsigned int > fNullList
Definition: dEdx_module.cc:89
QTextStream & endl(QTextStream &s)
Event finding and building.
std::vector< double > fMCdEdxVec
Definition: dEdx_module.cc:75
Signal from collection planes.
Definition: geo_types.h:146
void AnalysisExample::dEdx::beginJob ( )
virtual

Reimplemented from art::EDAnalyzer.

Definition at line 126 of file dEdx_module.cc.

127  {
128 
129  // hard-coded fd4apa cryostat dimensions
130  fStartEdepHist = tfs->make<TH3D>("fStartEdepHist","fStartEdepHist",50,-240.,240.,50,-750.,800.,50,-50.,560.);
131  fTrackEdepHist = tfs->make<TH3D>("fTrackEdepHist","fTrackEdepHist",50,-240.,240.,50,-750.,800.,50,-50.,560.);
132  fShowerEdepHist = tfs->make<TH3D>("fShowerEdepHist","fShowerEdepHist",50,-240.,240.,50,-750.,800.,50,-50.,560.);
133 
134  fStartEdepHist->GetXaxis()->SetTitle("x");
135  fStartEdepHist->GetYaxis()->SetTitle("y");
136  fStartEdepHist->GetZaxis()->SetTitle("z");
137 
138  fShowerEdepHist->GetXaxis()->SetTitle("x");
139  fShowerEdepHist->GetYaxis()->SetTitle("y");
140  fShowerEdepHist->GetZaxis()->SetTitle("z");
141 
142  fTrackEdepHist->GetXaxis()->SetTitle("x");
143  fTrackEdepHist->GetYaxis()->SetTitle("y");
144  fTrackEdepHist->GetZaxis()->SetTitle("z");
145 
146  fConvDistHist = tfs->make<TH1D>("fConvDistHist","fConvDistHist",100,0.,50.);
147 
148  }
art::ServiceHandle< art::TFileService > tfs
Definition: dEdx_module.cc:69
void AnalysisExample::dEdx::beginRun ( const art::Run run)

Definition at line 152 of file dEdx_module.cc.

153  {
154  }
void AnalysisExample::dEdx::checkDepositionVertex ( std::vector< sim::SimChannel SCHandle,
double  x,
double  y,
double  z 
)

Definition at line 275 of file dEdx_module.cc.

276  {
277  std::cout << "Checking for Deposition Vertex" << std::endl;
278  nearDeps = 0;
279  double track = 0;
280  for(auto const& channel : SCHandle)
281  {
282  if(fGeom->SignalType(channel.Channel()) == geo::kCollection)
283  {
284  auto const& timeSlices = channel.TDCIDEMap();
285  for(auto const& t : timeSlices)
286  {
287  auto const& eDeps = t.second;
288  for(auto const& eDep : eDeps)
289  {
290  track = std::sqrt( (eDep.x-x)*(eDep.x-x) + (eDep.y-y)*(eDep.y-y) + (eDep.z-z)*(eDep.z-z) );
291  if(track < 2.5 && eDep.energy > 0.)
292  {
293  nearDeps++;
294  } // if z position is minimum
295  } // energy deposition loop
296  } // time slice loop
297  } // if collection wire
298  } // sim channel loop
299  if( nearDeps > 30) goodDepVertex = true;
300  else ziFalse = z + 0.05;
301  } // checkDepositionVertex
art::ServiceHandle< geo::Geometry > fGeom
Definition: dEdx_module.cc:70
unsigned int nearDeps
Definition: dEdx_module.cc:98
uint8_t channel
Definition: CRTFragment.hh:201
SigType_t SignalType(geo::PlaneID const &pid) const
Returns the type of signal on the channels of specified TPC plane.
list x
Definition: train.py:276
QTextStream & endl(QTextStream &s)
Signal from collection planes.
Definition: geo_types.h:146
void AnalysisExample::dEdx::endJob ( )
virtual

Reimplemented from art::EDAnalyzer.

Definition at line 235 of file dEdx_module.cc.

236  {
237  fMCdEdxHist = tfs->make<TH1D>("fMCdEdxHist",";dE/dx (MeV/cm)",200,0.0,10.);
238 
239  for(unsigned int e = 0; e < fMCdEdxVec.size(); e++)
240  {
241  fMCdEdxHist->Fill(fMCdEdxVec[e]);
242  }
243 
244  std::cout << "Events with < 0.5 dE/dx entry: ";
245  for(unsigned int a = 0; a < fNullList.size(); a++) std::cout << fNullList[a] << ", ";
246 
247  }
const double e
art::ServiceHandle< art::TFileService > tfs
Definition: dEdx_module.cc:69
const double a
std::vector< unsigned int > fNullList
Definition: dEdx_module.cc:89
std::vector< double > fMCdEdxVec
Definition: dEdx_module.cc:75
void AnalysisExample::dEdx::findDepositionVertex ( std::vector< sim::SimChannel SCHandle)

Definition at line 249 of file dEdx_module.cc.

250  {
251  zi = 1000.;
252  std::cout << "Finding Deposition Vertex... " << std::endl;
253  for(auto const& channel : SCHandle)
254  {
255  if(fGeom->SignalType(channel.Channel()) == geo::kCollection)
256  {
257  auto const& timeSlices = channel.TDCIDEMap();
258  for(auto const& t : timeSlices)
259  {
260  auto const& eDeps = t.second;
261  for(auto const& eDep : eDeps)
262  {
263  if( (eDep.z < zi) && (eDep.energy > 0.) && (eDep.z > ziFalse) )
264  {
265  xi = eDep.x;
266  yi = eDep.y;
267  zi = eDep.z;
268  } // if z position is a "valid" minimum
269  } // energy deposition loop
270  } // time slice loop
271  } // if collection wire
272  } // sim channel loop
273  } // findDepositionVertex
art::ServiceHandle< geo::Geometry > fGeom
Definition: dEdx_module.cc:70
uint8_t channel
Definition: CRTFragment.hh:201
SigType_t SignalType(geo::PlaneID const &pid) const
Returns the type of signal on the channels of specified TPC plane.
QTextStream & endl(QTextStream &s)
Signal from collection planes.
Definition: geo_types.h:146
void AnalysisExample::dEdx::reconfigure ( fhicl::ParameterSet const &  pset)

Definition at line 158 of file dEdx_module.cc.

159  {
160  fSimulationProducerLabel = p.get<std::string>("SimulationLabel");
161  }
std::string string
Definition: nybbler.cc:12
std::string fSimulationProducerLabel
Definition: dEdx_module.cc:71
p
Definition: test.py:223

Member Data Documentation

double AnalysisExample::dEdx::dE
private

Definition at line 77 of file dEdx_module.cc.

double AnalysisExample::dEdx::fConvDist
private

Definition at line 102 of file dEdx_module.cc.

TH1D* AnalysisExample::dEdx::fConvDistHist
private

Definition at line 103 of file dEdx_module.cc.

art::ServiceHandle<geo::Geometry> AnalysisExample::dEdx::fGeom
private

Definition at line 70 of file dEdx_module.cc.

double AnalysisExample::dEdx::fMaxMCdEdx = 0.
private

Definition at line 76 of file dEdx_module.cc.

TH1D* AnalysisExample::dEdx::fMCdEdxHist
private

Definition at line 74 of file dEdx_module.cc.

std::vector<double> AnalysisExample::dEdx::fMCdEdxVec
private

Definition at line 75 of file dEdx_module.cc.

std::vector<unsigned int> AnalysisExample::dEdx::fNullList
private

Definition at line 89 of file dEdx_module.cc.

TH3D* AnalysisExample::dEdx::fShowerEdepHist
private

Definition at line 86 of file dEdx_module.cc.

std::string AnalysisExample::dEdx::fSimulationProducerLabel
private

Definition at line 71 of file dEdx_module.cc.

TH3D* AnalysisExample::dEdx::fStartEdepHist
private

Definition at line 84 of file dEdx_module.cc.

TH3D* AnalysisExample::dEdx::fTrackEdepHist
private

Definition at line 85 of file dEdx_module.cc.

bool AnalysisExample::dEdx::goodDepVertex = false
private

Definition at line 99 of file dEdx_module.cc.

unsigned int AnalysisExample::dEdx::nearDeps = 0
private

Definition at line 98 of file dEdx_module.cc.

art::ServiceHandle<art::TFileService> AnalysisExample::dEdx::tfs
private

Definition at line 69 of file dEdx_module.cc.

double AnalysisExample::dEdx::trackLength
private

Definition at line 80 of file dEdx_module.cc.

double AnalysisExample::dEdx::Vxi
private

Definition at line 104 of file dEdx_module.cc.

double AnalysisExample::dEdx::Vyi
private

Definition at line 105 of file dEdx_module.cc.

double AnalysisExample::dEdx::Vzi
private

Definition at line 106 of file dEdx_module.cc.

double AnalysisExample::dEdx::xi = 0.
private

Definition at line 92 of file dEdx_module.cc.

double AnalysisExample::dEdx::yi = 0.
private

Definition at line 93 of file dEdx_module.cc.

double AnalysisExample::dEdx::zi = 1000.
private

Definition at line 94 of file dEdx_module.cc.

double AnalysisExample::dEdx::ziFalse = -500.
private

Definition at line 97 of file dEdx_module.cc.


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