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

Public Member Functions

 ECalibration (fhicl::ParameterSet const &p)
 
 ECalibration (ECalibration const &)=delete
 
 ECalibration (ECalibration &&)=delete
 
ECalibrationoperator= (ECalibration const &)=delete
 
ECalibrationoperator= (ECalibration &&)=delete
 
void analyze (art::Event const &e) override
 
void beginJob () override
 
void beginRun (const art::Run &run) override
 
void reconfigure (fhicl::ParameterSet const &p)
 
- 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 Member Functions

void ResetVars ()
 
bool Has (std::vector< int > v, int i) const
 
double GetEdepMC (art::Event const &e) const
 
double GetEdepHits (detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp) const
 
bool EvMCselect (art::Event const &e)
 
void Make_dEdx (std::vector< double > &dEdx, std::vector< double > &range, const std::vector< proto::bHitInfo > &hits, double rmax) const
 
int GetClosestBeamTrkId () const
 
void SetBeamWindowEnd ()
 

Private Attributes

geo::GeometryCore const * fGeometry
 
double fElectronsToGeV
 
TVector3 fBeamPos
 
double fX0
 
double fX1
 
double fY0
 
double fY1
 
double fZ0
 
double fZ1
 
double fThXZ
 
double fThYZ
 
std::vector< int > fPdg
 
int fSimPdg
 
int fFlip
 
int fBestview
 
bool fStopping
 
bool fDecaying
 
double fMaxRange
 
double fT0
 
TH1D * fMomentumStartHist
 
TTree * fTree
 
TTree * fDataTree
 
int fRun
 
int fEvent
 
double fEdep
 
double fEdepMC
 
double fEkin
 
int fTrkIdx
 
double fdEdx
 
double fRange
 
std::map< size_t, std::vector< proto::bHitInfo >[3] > fTrk2InfoMap
 
art::Handle< std::vector< recob::Track > > fTrkListHandle
 
std::vector< art::Ptr< simb::MCParticle > > fSimlist
 
std::vector< art::Ptr< recob::Hit > > fHitlist
 
std::vector< art::Ptr< recob::Track > > fTracklist
 
std::vector< art::Ptr< recob::Vertex > > fVertexlist
 
std::vector< art::Ptr< recob::Shower > > fShslist
 
std::string fSimulationLabel
 
std::string fHitsModuleLabel
 
std::string fClusterModuleLabel
 
std::string fTrackModuleLabel
 
std::string fShowerModuleLabel
 
std::string fVertexModuleLabel
 
std::string fCalorimetryModuleLabel
 
calo::CalorimetryAlg fCalorimetryAlg
 

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 77 of file ECalibration_module.cc.

Constructor & Destructor Documentation

proto::ECalibration::ECalibration ( fhicl::ParameterSet const &  p)
explicit

Definition at line 178 of file ECalibration_module.cc.

179  :
180  EDAnalyzer(p),
181  fCalorimetryAlg(p.get<fhicl::ParameterSet>("CalorimetryAlg"))
182  // More initializers here.
183 {
184  // get a pointer to the geometry service provider
186 
187  reconfigure(p);
188 }
calo::CalorimetryAlg fCalorimetryAlg
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.h:25
p
Definition: test.py:223
geo::GeometryCore const * fGeometry
void reconfigure(fhicl::ParameterSet const &p)
proto::ECalibration::ECalibration ( ECalibration const &  )
delete
proto::ECalibration::ECalibration ( ECalibration &&  )
delete

Member Function Documentation

void proto::ECalibration::analyze ( art::Event const &  e)
overridevirtual

Implements art::EDAnalyzer.

Definition at line 335 of file ECalibration_module.cc.

336 {
337 
338  // Implementation of required member function here.
339  ResetVars();
341 
342  fRun = e.run();
343  fEvent = e.id().event();
344 
345  if (!EvMCselect(e)) return;
346 
347  // reco
348  // hits
349  auto hitListHandle = e.getHandle< std::vector<recob::Hit> >(fHitsModuleLabel);
350  if (hitListHandle) art::fill_ptr_vector(fHitlist, hitListHandle);
351 
352  auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService>()->DataFor(e);
353  auto const detProp = art::ServiceHandle<detinfo::DetectorPropertiesService>()->DataFor(e, clockData);
354 
355  if (fHitlist.size())
356  {
357  fEdep = GetEdepHits(clockData, detProp);
358  fEdepMC = GetEdepMC(e);
359  }
360 
361  // vertices
362  auto vtxListHandle = e.getHandle< std::vector<recob::Vertex> >(fVertexModuleLabel);
363  if (vtxListHandle) art::fill_ptr_vector(fVertexlist, vtxListHandle);
364 
365  // showers
366  auto shsListHandle = e.getHandle< std::vector<recob::Shower> >(fShowerModuleLabel);
367  if (shsListHandle) art::fill_ptr_vector(fShslist, shsListHandle);
368 
369  // tracks
370  fTrk2InfoMap.clear();
371  fTrkListHandle = e.getHandle< std::vector<recob::Track> >(fTrackModuleLabel);
372  if (fTrkListHandle)
373  {
374  art::FindManyP< recob::Hit, recob::TrackHitMeta > hitFromTrk(fTrkListHandle, e, fTrackModuleLabel);
375 
376  // choose the primary track using information about the beam window position
377  int idt = GetClosestBeamTrkId();
378  if (hitFromTrk.size() && (idt > -1))
379  {
380  fFlip = 0; // tag track if it is reconstructed in the opposite direction
381  if ((*fTrkListHandle)[idt].End().Z() < (*fTrkListHandle)[idt].Vertex().Z())
382  { fFlip = 1; }
383 
384  auto vhit = hitFromTrk.at(idt);
385  auto vmeta = hitFromTrk.data(idt);
386 
387  for (size_t h = 0; h < vhit.size(); ++h)
388  {
389  int view = vhit[h]->WireID().Plane;
390 
391  if (view != fBestview) continue;
392 
393  size_t idx = vmeta[h]->Index();
394  double tdrift = vhit[h]->PeakTime();
395  double dx = vmeta[h]->Dx();
396  double dqadc = vhit[h]->Integral();
397  int wire = vhit[h]->WireID().Wire;
398 
399  double dq = fCalorimetryAlg.dEdx_AREA(clockData, detProp, dqadc/dx, tdrift, view, fT0) * dx;
400  fEkin += dq;
401 
402  fTrk2InfoMap[idt][view].emplace_back(idx, dx, dq, wire);
403  }
404 
405  for (auto const & trkEntry : fTrk2InfoMap)
406  {
407  std::vector< double > dEdx, range;
408  auto const & info = trkEntry.second;
409  Make_dEdx(dEdx, range, info[fBestview], fMaxRange);
410 
411  for (size_t i = 0; i < dEdx.size(); ++i)
412  {
413  fTrkIdx = trkEntry.first;
414  fdEdx = dEdx[i];
415  fRange = range[i];
416 
417  fDataTree->Fill();
418  }
419  }
420  }
421  }
422  fTree->Fill();
423 }
calo::CalorimetryAlg fCalorimetryAlg
std::vector< art::Ptr< recob::Hit > > fHitlist
art::Handle< std::vector< recob::Track > > fTrkListHandle
int GetClosestBeamTrkId() const
double dEdx(float dqdx, float Efield)
Definition: doAna.cpp:21
const double e
double GetEdepHits(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp) const
std::vector< art::Ptr< recob::Vertex > > fVertexlist
void Make_dEdx(std::vector< double > &dEdx, std::vector< double > &range, const std::vector< proto::bHitInfo > &hits, double rmax) const
double GetEdepMC(art::Event const &e) const
double dEdx_AREA(detinfo::DetectorClocksData const &clock_data, detinfo::DetectorPropertiesData const &det_prop, recob::Hit const &hit, double pitch, double T0=0) const
bool EvMCselect(art::Event const &e)
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:297
std::vector< art::Ptr< recob::Shower > > fShslist
std::map< size_t, std::vector< proto::bHitInfo >[3] > fTrk2InfoMap
void proto::ECalibration::beginJob ( )
overridevirtual

Reimplemented from art::EDAnalyzer.

Definition at line 196 of file ECalibration_module.cc.

197 {
198  // access art's TFileService, which will handle creating and writing hists
200 
201  fMomentumStartHist = tfs->make<TH1D>("mom",";particle Momentum (GeV);",100, 0., 0.);
202 
203  fTree = tfs->make<TTree>("calibration","calibration tree");
204  fTree->Branch("fRun", &fRun, "fRun/I");
205  fTree->Branch("fEvent", &fEvent, "fEvent/I");
206  fTree->Branch("fEdep", &fEdep, "fEdep/D");
207  fTree->Branch("fEdepMC", &fEdepMC, "fEdepMC/D");
208  fTree->Branch("fEkin", &fEkin, "fEkin/D");
209 
210  fDataTree = tfs->make<TTree>("Data", "dE/dx info");
211  fDataTree->Branch("fRun", &fRun, "fRun/I");
212  fDataTree->Branch("fEvent", &fEvent, "fEvent/I");
213  fDataTree->Branch("fTrkIdx", &fTrkIdx, "fTrkIdx/I");
214  fDataTree->Branch("fdEdx", &fdEdx, "fdEdx/D");
215  fDataTree->Branch("fRange", &fRange, "fRange/D");
216 }
void proto::ECalibration::beginRun ( const art::Run run)
override

Definition at line 190 of file ECalibration_module.cc.

191 {
193  fElectronsToGeV = 1./larParameters->GeVToElectrons();
194 }
double GeVToElectrons() const
bool proto::ECalibration::EvMCselect ( art::Event const &  e)
private

Definition at line 246 of file ECalibration_module.cc.

247 {
248  std::vector< art::Ptr<simb::MCParticle> > simlist;
249 
250  auto mcparticleHandle = e.getHandle< std::vector<simb::MCParticle> >(fSimulationLabel);
251  if (mcparticleHandle)
252  art::fill_ptr_vector(simlist, mcparticleHandle);
253 
254  std::map< int, const simb::MCParticle* > particleMap;
255  for (auto const& particle : simlist)
256  {
257  particleMap[particle->TrackId()] = &*particle;
258  }
259 
260  fSimPdg = particleMap.begin()->second->PdgCode();
261  fT0 = particleMap.begin()->second->T();
262 
263  if (fStopping)
264  {
265  if ((fSimPdg == 2212)
266  && Has(fPdg, fSimPdg)
267  && (particleMap.size() == 1)
268  && (particleMap.begin()->second->EndProcess() != "ProtonInelastic"))
269  {
270  return true;
271  }
272  else if ((fSimPdg == 13)
273  && Has(fPdg, fSimPdg))
274  {
275  return true;
276  }
277  else if ((fSimPdg == 211)
278  && (Has(fPdg, fSimPdg)))
279  {
280  return true;
281  }
282  }
283 
284  if (fDecaying)
285  {
286  if ((fSimPdg == 321)
287  && Has(fPdg, fSimPdg)
288  && (particleMap.begin()->second->EndProcess() == "Decay")
289  && (particleMap.begin()->second->NumberDaughters() == 2)
290  && (particleMap.size() == 6))
291  {
292  bool kaonmodedcy = false;
293  size_t count = 0;
294  for (auto const & p : particleMap)
295  {
296  if ((p.second->PdgCode() == 321) ||
297  (p.second->PdgCode() == -13) ||
298  (p.second->PdgCode() == 14))
299  {
300  kaonmodedcy = true;
301  }
302 
303  count++;
304  if (count == 3) break;
305  }
306  return kaonmodedcy;
307  }
308  else if ((fSimPdg == 13)
309  && Has(fPdg, fSimPdg))
310  {
311  return true;
312  }
313  else if ((fSimPdg == 211)
314  && (Has(fPdg, fSimPdg))
315  && (particleMap.begin()->second->EndProcess() == "Decay"))
316  {
317  return true;
318  }
319  }
320 
321  if (!fDecaying
322  && !fStopping
323  && (fSimPdg == abs(11))
324  && Has(fPdg, fSimPdg))
325  {
326  return true;
327  }
328 
329  fMomentumStartHist->Fill(particleMap.begin()->second->Momentum(0).P());
330 
331  return false;
332 }
T abs(T value)
const double e
p
Definition: test.py:223
std::vector< int > fPdg
bool Has(std::vector< int > v, int i) const
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:297
int proto::ECalibration::GetClosestBeamTrkId ( ) const
private

Definition at line 425 of file ECalibration_module.cc.

426 {
427  int idt = -1;
428  if (!fTrkListHandle->size()) {return idt;}
429 
430  auto const & trk0 = (*fTrkListHandle)[0];
431  float mindist2 = pma::Dist2(fBeamPos, trk0.Vertex());
432  if (trk0.End().Z() < trk0.Vertex().Z())
433  {
434  mindist2 = pma::Dist2(fBeamPos, trk0.End());
435  idt = 0;
436  }
437  else
438  {
439  mindist2 = pma::Dist2(fBeamPos, trk0.Vertex());
440  idt = 0;
441  }
442 
443  for (size_t t = 1; t < fTrkListHandle->size(); ++t)
444  {
445  auto const & trk = (*fTrkListHandle)[t];
446 
447  float dist2 = pma::Dist2(fBeamPos, trk.Vertex());
448  if (trk.End().Z() < trk.Vertex().Z())
449  {
450  dist2 = pma::Dist2(fBeamPos, trk.End());
451  }
452  else
453  {
454  dist2 = pma::Dist2(fBeamPos, trk.Vertex());
455  }
456 
457  //
458  // idt:
459  // we assumed that primary, incoming particle
460  // has an index of a track closest to the beam window
461  //
462  if (dist2 < mindist2) {mindist2 = dist2; idt = t;}
463  }
464 
465  return idt;
466 }
art::Handle< std::vector< recob::Track > > fTrkListHandle
double Dist2(const TVector2 &v1, const TVector2 &v2)
Definition: Utilities.cxx:37
double proto::ECalibration::GetEdepHits ( detinfo::DetectorClocksData const &  clockData,
detinfo::DetectorPropertiesData const &  detProp 
) const
private

Definition at line 533 of file ECalibration_module.cc.

535 {
536  if (!fHitlist.size()) return 0.0;
537 
538  double dqsum = 0.0;
539  for (size_t h = 0; h < fHitlist.size(); ++h)
540  {
541  unsigned short plane = fHitlist[h]->WireID().Plane;
542  if (plane != geo::kZ) continue;
543 
544  double dqadc = fHitlist[h]->Integral();
545  if (!std::isnormal(dqadc) || (dqadc < 0)) continue;
546 
547  double dqel = fCalorimetryAlg.ElectronsFromADCArea(dqadc, plane);
548 
549  double tdrift = fHitlist[h]->PeakTime();
550  double correllifetime = fCalorimetryAlg.LifetimeCorrection(clockData, detProp, tdrift, fT0);
551 
552  double dq = dqel * correllifetime * fElectronsToGeV * 1000;
553  if (!std::isnormal(dq) || (dq < 0)) continue;
554 
555  dqsum += dq;
556  }
557 
558  return dqsum;
559 }
calo::CalorimetryAlg fCalorimetryAlg
std::vector< art::Ptr< recob::Hit > > fHitlist
Planes which measure Z direction.
Definition: geo_types.h:132
double ElectronsFromADCArea(double area, unsigned short plane) const
double LifetimeCorrection(detinfo::DetectorClocksData const &clock_data, detinfo::DetectorPropertiesData const &det_prop, double time, double T0=0) const
double proto::ECalibration::GetEdepMC ( art::Event const &  e) const
private

Definition at line 505 of file ECalibration_module.cc.

506 {
507  double energy = 0.0;
508 
509  auto simchannelHandle = e.getHandle< std::vector<sim::SimChannel> >(fSimulationLabel);
510  if (simchannelHandle)
511  {
512  for ( auto const& channel : (*simchannelHandle) )
513  {
514  // for every time slice in this channel:
515  auto const& timeSlices = channel.TDCIDEMap();
516  for ( auto const& timeSlice : timeSlices )
517  {
518  // loop over the energy deposits.
519  auto const& energyDeposits = timeSlice.second;
520 
521  for ( auto const& energyDeposit : energyDeposits )
522  {
523  energy += energyDeposit.numElectrons * fElectronsToGeV * 1000;
524  }
525  }
526  }
527  }
528 
529  return energy;
530 }
uint8_t channel
Definition: CRTFragment.hh:201
const double e
bool proto::ECalibration::Has ( std::vector< int >  v,
int  i 
) const
private

Definition at line 561 of file ECalibration_module.cc.

562 {
563  for (auto c : v) if (c == i) return true;
564  return false;
565 }
void proto::ECalibration::Make_dEdx ( std::vector< double > &  dEdx,
std::vector< double > &  range,
const std::vector< proto::bHitInfo > &  hits,
double  rmax 
) const
private

Definition at line 468 of file ECalibration_module.cc.

469 {
470  if (!hits.size()) return;
471 
472  dEdx.clear(); range.clear();
473 
474  int i0 = hits.size() - 1; int i1 = -1; int di = -1;
475  if (Has(fPdg, abs(11)) || fFlip) {i0 = 0; i1 = hits.size(); di = 1;}
476 
477  double de = 0.0;
478  double dx = 0.0;
479  double r0 = 0.0; double r1 = 0.0; double r = 0.0;
480 
481  double minDx = 0.1; // can be a parameter
482 
483  while ((i0 != i1) && (r < rmax))
484  {
485  dx = 0.0; de = 0.0;
486  while ((i0 != i1) && (dx <= minDx))
487  {
488  de += hits[i0].dE;
489  dx += hits[i0].dx;
490  i0 += di;
491  }
492 
493  r0 = r1;
494  r1 += dx;
495  r = 0.5 * (r0 + r1);
496 
497  if ((de > 0.0) && (dx > 0.0) && (r < rmax))
498  {
499  dEdx.push_back(de/dx);
500  range.push_back(r);
501  }
502  }
503 }
T abs(T value)
std::vector< int > fPdg
bool Has(std::vector< int > v, int i) const
ECalibration& proto::ECalibration::operator= ( ECalibration const &  )
delete
ECalibration& proto::ECalibration::operator= ( ECalibration &&  )
delete
void proto::ECalibration::reconfigure ( fhicl::ParameterSet const &  p)

Definition at line 218 of file ECalibration_module.cc.

219 {
220  fSimulationLabel = p.get< std::string >("SimulationLabel");
221  fHitsModuleLabel = p.get< std::string >("HitsModuleLabel");
222  fTrackModuleLabel = p.get< std::string >("TrackModuleLabel");
223  fShowerModuleLabel = p.get< std::string >("ShowerModuleLabel");
224  fClusterModuleLabel = p.get< std::string >("ClusterModuleLabel");
225  fVertexModuleLabel = p.get< std::string >("VertexModuleLabel");
226  fCalorimetryModuleLabel = p.get< std::string >("CalorimetryModuleLabel");
227  fX0 = p.get<double>("BeamPosX");
228  fY0 = p.get<double>("BeamPosY");
229  fZ0 = p.get<double>("BeamPosZ");
230  fThXZ = p.get<double>("ThXZ");
231  fThYZ = p.get<double>("ThYZ");
232 
233  fStopping = p.get<bool>("Stopping");
234  fDecaying = p.get<bool>("Decaying");
235 
236  fMaxRange = p.get<double>("MaxRange");
237  fPdg = p.get< std::vector<int> >("Pdg");
238 
239  fBestview = p.get<int>("Bestview");
240 }
std::string string
Definition: nybbler.cc:12
p
Definition: test.py:223
std::vector< int > fPdg
std::string fCalorimetryModuleLabel
void proto::ECalibration::ResetVars ( void  )
private

Definition at line 582 of file ECalibration_module.cc.

583 {
584  fSimlist.clear();
585  fHitlist.clear();
586  fTracklist.clear();
587  fVertexlist.clear();
588  fShslist.clear();
589  fSimPdg = 0;
590  fEdep = 0.0;
591  fEdepMC = 0.0;
592  fEkin = 0;
593  fT0 = 0.0;
594  fTrkIdx = 0;
595  fRange = 0.0;
596  fdEdx = 0.0;
597  fFlip = 0;
598 }
std::vector< art::Ptr< recob::Hit > > fHitlist
std::vector< art::Ptr< simb::MCParticle > > fSimlist
std::vector< art::Ptr< recob::Vertex > > fVertexlist
std::vector< art::Ptr< recob::Track > > fTracklist
std::vector< art::Ptr< recob::Shower > > fShslist
void proto::ECalibration::SetBeamWindowEnd ( )
private

Definition at line 568 of file ECalibration_module.cc.

569 {
570  fZ1 = 0.0;
571 
572  double dz = fZ1 - fZ0;
573  double dx = dz * tan(fThXZ * (TMath::Pi() / 180.0));
574  double dy = dz * tan(fThYZ * (TMath::Pi() / 180.0));
575 
576  fX1 = fX0 + dx;
577  fY1 = fY0 + dy;
578 
579  fBeamPos.SetXYZ(fX1, fY1, fZ1);
580 }

Member Data Documentation

TVector3 proto::ECalibration::fBeamPos
private

Definition at line 122 of file ECalibration_module.cc.

int proto::ECalibration::fBestview
private

Definition at line 134 of file ECalibration_module.cc.

calo::CalorimetryAlg proto::ECalibration::fCalorimetryAlg
private

Definition at line 174 of file ECalibration_module.cc.

std::string proto::ECalibration::fCalorimetryModuleLabel
private

Definition at line 172 of file ECalibration_module.cc.

std::string proto::ECalibration::fClusterModuleLabel
private

Definition at line 168 of file ECalibration_module.cc.

TTree* proto::ECalibration::fDataTree
private

Definition at line 145 of file ECalibration_module.cc.

bool proto::ECalibration::fDecaying
private

Definition at line 137 of file ECalibration_module.cc.

double proto::ECalibration::fdEdx
private

Definition at line 152 of file ECalibration_module.cc.

double proto::ECalibration::fEdep
private

Definition at line 148 of file ECalibration_module.cc.

double proto::ECalibration::fEdepMC
private

Definition at line 149 of file ECalibration_module.cc.

double proto::ECalibration::fEkin
private

Definition at line 150 of file ECalibration_module.cc.

double proto::ECalibration::fElectronsToGeV
private

Definition at line 103 of file ECalibration_module.cc.

int proto::ECalibration::fEvent
private

Definition at line 147 of file ECalibration_module.cc.

int proto::ECalibration::fFlip
private

Definition at line 133 of file ECalibration_module.cc.

geo::GeometryCore const* proto::ECalibration::fGeometry
private

Definition at line 102 of file ECalibration_module.cc.

std::vector< art::Ptr<recob::Hit> > proto::ECalibration::fHitlist
private

Definition at line 160 of file ECalibration_module.cc.

std::string proto::ECalibration::fHitsModuleLabel
private

Definition at line 167 of file ECalibration_module.cc.

double proto::ECalibration::fMaxRange
private

Definition at line 139 of file ECalibration_module.cc.

TH1D* proto::ECalibration::fMomentumStartHist
private

Definition at line 143 of file ECalibration_module.cc.

std::vector<int> proto::ECalibration::fPdg
private

Definition at line 130 of file ECalibration_module.cc.

double proto::ECalibration::fRange
private

Definition at line 153 of file ECalibration_module.cc.

int proto::ECalibration::fRun
private

Definition at line 146 of file ECalibration_module.cc.

std::string proto::ECalibration::fShowerModuleLabel
private

Definition at line 170 of file ECalibration_module.cc.

std::vector< art::Ptr<recob::Shower> > proto::ECalibration::fShslist
private

Definition at line 163 of file ECalibration_module.cc.

std::vector< art::Ptr<simb::MCParticle> > proto::ECalibration::fSimlist
private

Definition at line 159 of file ECalibration_module.cc.

int proto::ECalibration::fSimPdg
private

Definition at line 131 of file ECalibration_module.cc.

std::string proto::ECalibration::fSimulationLabel
private

Definition at line 166 of file ECalibration_module.cc.

bool proto::ECalibration::fStopping
private

Definition at line 136 of file ECalibration_module.cc.

double proto::ECalibration::fT0
private

Definition at line 140 of file ECalibration_module.cc.

double proto::ECalibration::fThXZ
private

Definition at line 126 of file ECalibration_module.cc.

double proto::ECalibration::fThYZ
private

Definition at line 127 of file ECalibration_module.cc.

std::vector< art::Ptr<recob::Track> > proto::ECalibration::fTracklist
private

Definition at line 161 of file ECalibration_module.cc.

std::string proto::ECalibration::fTrackModuleLabel
private

Definition at line 169 of file ECalibration_module.cc.

TTree* proto::ECalibration::fTree
private

Definition at line 145 of file ECalibration_module.cc.

std::map< size_t, std::vector< proto::bHitInfo >[3] > proto::ECalibration::fTrk2InfoMap
private

Definition at line 156 of file ECalibration_module.cc.

int proto::ECalibration::fTrkIdx
private

Definition at line 151 of file ECalibration_module.cc.

art::Handle< std::vector<recob::Track> > proto::ECalibration::fTrkListHandle
private

Definition at line 157 of file ECalibration_module.cc.

std::vector< art::Ptr<recob::Vertex> > proto::ECalibration::fVertexlist
private

Definition at line 162 of file ECalibration_module.cc.

std::string proto::ECalibration::fVertexModuleLabel
private

Definition at line 171 of file ECalibration_module.cc.

double proto::ECalibration::fX0
private

Definition at line 123 of file ECalibration_module.cc.

double proto::ECalibration::fX1
private

Definition at line 123 of file ECalibration_module.cc.

double proto::ECalibration::fY0
private

Definition at line 124 of file ECalibration_module.cc.

double proto::ECalibration::fY1
private

Definition at line 124 of file ECalibration_module.cc.

double proto::ECalibration::fZ0
private

Definition at line 125 of file ECalibration_module.cc.

double proto::ECalibration::fZ1
private

Definition at line 125 of file ECalibration_module.cc.


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