Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
10 changes: 4 additions & 6 deletions sbncode/OpT0Finder/flashmatch/Algorithms/ChargeAnalytical.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,6 @@

#include "ChargeAnalytical.h"

#pragma GCC diagnostic ignored "-Wdangling-reference"

namespace flashmatch {

static ChargeAnalyticalFactory __global_ChargeAnalyticalFactory__;
Expand All @@ -17,9 +15,9 @@ namespace flashmatch {
{
_global_qe = pset.get<double>("GlobalQE");
_qe_v = pset.get<std::vector<double> >("VUVEfficiency");
if(_qe_v.size() != DetectorSpecs::GetME().NOpDets()) {
if(_qe_v.size() != DetectorSpecs::GetME()->NOpDets()) {
FLASH_CRITICAL() << "VUV Efficiency array has size " << _qe_v.size()
<< " != number of opdet (" << DetectorSpecs::GetME().NOpDets() << ")!" << std::endl;
<< " != number of opdet (" << DetectorSpecs::GetME()->NOpDets() << ")!" << std::endl;
throw OpT0FinderException();
}
}
Expand All @@ -28,7 +26,7 @@ namespace flashmatch {
Flash_t &flash) const
{

size_t n_pmt = DetectorSpecs::GetME().NOpDets();
size_t n_pmt = DetectorSpecs::GetME()->NOpDets();

for (size_t i = 0; i < n_pmt; ++i) {
flash.pe_v[i] = 0;
Expand All @@ -40,7 +38,7 @@ namespace flashmatch {

auto const &pt = track[pt_index];

auto const& pmt_pos = DetectorSpecs::GetME().PMTPosition(pmt_index);
auto const& pmt_pos = DetectorSpecs::GetME()->PMTPosition(pmt_index);
double dx = pmt_pos[0] - pt.x;
double dy = pmt_pos[1] - pt.y;
double dz = pmt_pos[2] - pt.z;
Expand Down
4 changes: 1 addition & 3 deletions sbncode/OpT0Finder/flashmatch/Algorithms/LightPath.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,6 @@

#include "LightPath.h"

#pragma GCC diagnostic ignored "-Wdangling-reference"

namespace flashmatch {

static LightPathFactory __global_LightPathFactory__;
Expand Down Expand Up @@ -87,7 +85,7 @@ namespace flashmatch {
}

// Trimming Q_cluster
auto const& bbox = DetectorSpecs::GetME().ActiveVolume();
auto const& bbox = DetectorSpecs::GetME()->ActiveVolume();
double _vol_xmax = bbox.Max()[0];
double _vol_ymax = bbox.Max()[1];
double _vol_zmax = bbox.Max()[2];
Expand Down
26 changes: 13 additions & 13 deletions sbncode/OpT0Finder/flashmatch/Algorithms/PhotonLibHypothesis.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -34,25 +34,25 @@ namespace flashmatch {
_qe_refl_v.clear();
_qe_v = pset.get<std::vector<double> >("VUVEfficiency",_qe_v);
_qe_refl_v = pset.get<std::vector<double> >("VISEfficiency",_qe_refl_v);
if(_qe_v.empty()) _qe_v.resize(DetectorSpecs::GetME().NOpDets(),1.0);
if(_qe_refl_v.empty()) _qe_refl_v.resize(DetectorSpecs::GetME().NOpDets(),1.0);
if(_qe_v.size() != DetectorSpecs::GetME().NOpDets()) {
if(_qe_v.empty()) _qe_v.resize(DetectorSpecs::GetME()->NOpDets(),1.0);
if(_qe_refl_v.empty()) _qe_refl_v.resize(DetectorSpecs::GetME()->NOpDets(),1.0);
if(_qe_v.size() != DetectorSpecs::GetME()->NOpDets()) {
FLASH_CRITICAL() << "VUV Efficiency factor array has size " << _qe_v.size()
<< " != number of opdet (" << DetectorSpecs::GetME().NOpDets() << ")!" << std::endl;
<< " != number of opdet (" << DetectorSpecs::GetME()->NOpDets() << ")!" << std::endl;
throw OpT0FinderException();
}
if(_qe_refl_v.size() != DetectorSpecs::GetME().NOpDets()) {
if(_qe_refl_v.size() != DetectorSpecs::GetME()->NOpDets()) {
FLASH_CRITICAL() << "VIS Efficiency factor array has size " << _qe_refl_v.size()
<< " != number of opdet (" << DetectorSpecs::GetME().NOpDets() << ")!" << std::endl;
<< " != number of opdet (" << DetectorSpecs::GetME()->NOpDets() << ")!" << std::endl;
throw OpT0FinderException();
}

// By default, add all opdets to the channel mask
// Note that this may be overridden by the manager
// via the SetChannelMask() method.
_channel_mask.clear();
_channel_mask.reserve(DetectorSpecs::GetME().NOpDets());
for (size_t i = 0; i < DetectorSpecs::GetME().NOpDets(); i++) {
_channel_mask.reserve(DetectorSpecs::GetME()->NOpDets());
for (size_t i = 0; i < DetectorSpecs::GetME()->NOpDets(); i++) {
_channel_mask[i] = i;
}
}
Expand Down Expand Up @@ -140,7 +140,7 @@ namespace flashmatch {

static double xyz[3] = {0.};

size_t n_pmt = DetectorSpecs::GetME().NOpDets();
size_t n_pmt = DetectorSpecs::GetME()->NOpDets();

for ( size_t ipmt = 0; ipmt < n_pmt; ++ipmt) {

Expand Down Expand Up @@ -193,12 +193,12 @@ namespace flashmatch {

static double xyz[3] = {0.};

size_t n_pmt = DetectorSpecs::GetME().NOpDets();
size_t n_pmt = DetectorSpecs::GetME()->NOpDets();

for ( auto& v : flash.pe_v ) v = 0;


size_t n_pmt = DetectorSpecs::GetME().NOpDets();//n_pmt returns 0 now, needs to be fixed
size_t n_pmt = DetectorSpecs::GetME()->NOpDets();//n_pmt returns 0 now, needs to be fixed
if(flash.pe_v.empty()) flash.pe_v.resize(n_pmt);
if(flash.pe_err_v.empty()) flash.pe_err_v.resize(n_pmt);

Expand All @@ -210,7 +210,7 @@ namespace flashmatch {

auto det = DetectorSpecs::GetME();

auto const& lib_data = DetectorSpecs::GetME().GetPhotonLibraryData();
auto const& lib_data = DetectorSpecs::GetME()->GetPhotonLibraryData();

//start = high_resolution_clock::now();
#pragma omp parallel
Expand All @@ -221,7 +221,7 @@ namespace flashmatch {
size_t start_pt = num_pts * thread_id;
if(thread_id+1 == num_threads) num_pts += (trk.size() % num_threads);

auto const& vox_def = DetectorSpecs::GetME().GetVoxelDef();
auto const& vox_def = DetectorSpecs::GetME()->GetVoxelDef();
// auto s = vox_def.GetVoxelSize();
// auto s1 = vox_def.GetRegionLowerCorner();
// auto s2 = vox_def.GetRegionUpperCorner();
Expand Down
48 changes: 23 additions & 25 deletions sbncode/OpT0Finder/flashmatch/Algorithms/QLLMatch.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,6 @@

#include "QLLMatch.h"

#pragma GCC diagnostic ignored "-Wdangling-reference"

using namespace std::chrono;
namespace flashmatch {

Expand Down Expand Up @@ -48,21 +46,21 @@ namespace flashmatch {
_onepmt_pesum_threshold = pset.get<double>("OnePMTPESumThreshold");
_onepmt_pefrac_threshold = pset.get<double>("OnePMTPEFracThreshold");

_xpos_v.resize(DetectorSpecs::GetME().NOpDets(),0.);
_ypos_v.resize(DetectorSpecs::GetME().NOpDets(),0.);
_zpos_v.resize(DetectorSpecs::GetME().NOpDets(),0.);
for(size_t ch=0; ch<DetectorSpecs::GetME().NOpDets(); ++ch) {
auto const& pmt_pos = DetectorSpecs::GetME().PMTPosition(ch);
_xpos_v.resize(DetectorSpecs::GetME()->NOpDets(),0.);
_ypos_v.resize(DetectorSpecs::GetME()->NOpDets(),0.);
_zpos_v.resize(DetectorSpecs::GetME()->NOpDets(),0.);
for(size_t ch=0; ch<DetectorSpecs::GetME()->NOpDets(); ++ch) {
auto const& pmt_pos = DetectorSpecs::GetME()->PMTPosition(ch);
_xpos_v[ch] = pmt_pos[0];
_ypos_v[ch] = pmt_pos[1];
_zpos_v[ch] = pmt_pos[2];
}
if(_tpc == -1 || _cryo == -1) {
auto const& bbox = DetectorSpecs::GetME().ActiveVolume();
auto const& bbox = DetectorSpecs::GetME()->ActiveVolume();
_vol_xmax = bbox.Max()[0];
_vol_xmin = bbox.Min()[0];
} else {
auto const& bbox = DetectorSpecs::GetME().ActiveVolume(_tpc, _cryo);
auto const& bbox = DetectorSpecs::GetME()->ActiveVolume(_tpc, _cryo);
_vol_xmax = bbox.Max()[0];
_vol_xmin = bbox.Min()[0];
}
Expand All @@ -71,8 +69,8 @@ namespace flashmatch {
// Note that this may be overridden by the manager
// via the SetChannelMask() method.
_channel_mask.clear();
_channel_mask.reserve(DetectorSpecs::GetME().NOpDets());
for (size_t i = 0; i < DetectorSpecs::GetME().NOpDets(); i++) {
_channel_mask.reserve(DetectorSpecs::GetME()->NOpDets());
for (size_t i = 0; i < DetectorSpecs::GetME()->NOpDets(); i++) {
_channel_mask[i] = i;
}
}
Expand All @@ -81,7 +79,7 @@ namespace flashmatch {
_tpc = tpc;
_cryo = cryo;

auto const& bbox = DetectorSpecs::GetME().ActiveVolume(_tpc, _cryo);
auto const& bbox = DetectorSpecs::GetME()->ActiveVolume(_tpc, _cryo);
_vol_xmax = bbox.Max()[0];
_vol_xmin = bbox.Min()[0];
}
Expand All @@ -93,7 +91,7 @@ namespace flashmatch {

// combine cluster + flash mask for this match pair
_match_mask.clear();
_match_mask.resize(DetectorSpecs::GetME().NOpDets(), 0);
_match_mask.resize(DetectorSpecs::GetME()->NOpDets(), 0);

for (size_t opch=0; opch < flash.pe_v.size(); opch++){
if (flash.pds_mask_v.at(opch)!=0 || pt_v.tpc_mask_v.at(opch)!=0){
Expand Down Expand Up @@ -150,15 +148,15 @@ namespace flashmatch {

// initialize the hypothesis flash
Flash_t one_hypothesis;
one_hypothesis.pe_v.resize(DetectorSpecs::GetME().NOpDets(), 0.);
one_hypothesis.pe_v.resize(DetectorSpecs::GetME()->NOpDets(), 0.);
for (auto &v : one_hypothesis.pe_v) v = 0;

FillEstimate(_raw_trk,one_hypothesis);

// initialize the measurement flash
auto one_measurement = flash;

for (size_t ich = 0; ich < DetectorSpecs::GetME().NOpDets(); ++ich ) {
for (size_t ich = 0; ich < DetectorSpecs::GetME()->NOpDets(); ++ich ) {
if (_match_mask.at(ich) != 0){
one_hypothesis.pe_v[ich] = 0.;
one_measurement.pe_v[ich] = 0.;
Expand All @@ -169,7 +167,7 @@ namespace flashmatch {
// - when the measured flash PE is equal to 0 and the hypothesis is large, assume that the
// measured flash PE was set to 0 due to saturation
if (_saturated_thresh > 0){
for (size_t ich = 0; ich < DetectorSpecs::GetME().NOpDets(); ++ich ) {
for (size_t ich = 0; ich < DetectorSpecs::GetME()->NOpDets(); ++ich ) {
// if above the saturated threshold, measured is zero, is a PMT, and is not masked
if ((one_hypothesis.pe_v[ich] >= _saturated_thresh) && (one_measurement.pe_v[ich] == 0) && (_channel_type[ich] == 0) && (_match_mask.at(ich) == 0)){
std::cout << "Guessing " << ich << " is saturated, setting hypothesis to 0" << std::endl;
Expand Down Expand Up @@ -235,7 +233,7 @@ FlashMatch_t QLLMatch::OnePMTMatch(const Flash_t& flash) {
}

// Now see if Flash T0 can be consistent with an assumption MinX @ X=0.
double xdiff = fabs(_raw_xmin_pt.x - flash.time * DetectorSpecs::GetME().DriftVelocity());
double xdiff = fabs(_raw_xmin_pt.x - flash.time * DetectorSpecs::GetME()->DriftVelocity());
if( xdiff > _onepmt_xdiff_threshold ) {
//std::cout << "XDiffThreshold not met (xdiff=" << xdiff << ")" << std::endl;
return res;
Expand All @@ -253,7 +251,7 @@ FlashMatch_t QLLMatch::OnePMTMatch(const Flash_t& flash) {
// Compute TPC point
res.tpc_point.x = res.tpc_point.y = res.tpc_point.z = 0;
double weight = 0;
for (size_t pmt_index = 0; pmt_index < DetectorSpecs::GetME().NOpDets(); ++pmt_index) {
for (size_t pmt_index = 0; pmt_index < DetectorSpecs::GetME()->NOpDets(); ++pmt_index) {

res.tpc_point.y += _ypos_v.at(pmt_index) * _hypothesis.pe_v[pmt_index];
res.tpc_point.z += _zpos_v.at(pmt_index) * _hypothesis.pe_v[pmt_index];
Expand Down Expand Up @@ -293,7 +291,7 @@ FlashMatch_t QLLMatch::OnePMTMatch(const Flash_t& flash) {

double weight = 0;

for (size_t pmt_index = 0; pmt_index < DetectorSpecs::GetME().NOpDets(); ++pmt_index) {
for (size_t pmt_index = 0; pmt_index < DetectorSpecs::GetME()->NOpDets(); ++pmt_index) {

res.tpc_point.y += _ypos_v.at(pmt_index) * _hypothesis.pe_v[pmt_index];
res.tpc_point.z += _zpos_v.at(pmt_index) * _hypothesis.pe_v[pmt_index];
Expand All @@ -319,13 +317,13 @@ FlashMatch_t QLLMatch::OnePMTMatch(const Flash_t& flash) {

// Compute X-weighting
/*
double x0 = _raw_xmin_pt.x - flash.time * DetectorSpecs::GetME().DriftVelocity();
double x0 = _raw_xmin_pt.x - flash.time * DetectorSpecs::GetME()->DriftVelocity();
if( fabs(_reco_x_offset - x0) > _recox_penalty_threshold )
res.score *= 1. / (1. + fabs(_reco_x_offset - x0) - _recox_penalty_threshold);
// Compute Z-weighting
double z0 = 0;
weight = 0;
for (size_t pmt_index = 0; pmt_index < DetectorSpecs::GetME().NOpDets(); ++pmt_index) {
for (size_t pmt_index = 0; pmt_index < DetectorSpecs::GetME()->NOpDets(); ++pmt_index) {
z0 += _zpos_v.at(pmt_index) * flash.pe_v[pmt_index];
weight += flash.pe_v[pmt_index];
}
Expand All @@ -338,8 +336,8 @@ FlashMatch_t QLLMatch::OnePMTMatch(const Flash_t& flash) {

const Flash_t &QLLMatch::ChargeHypothesis(const double xoffset) {
auto start = high_resolution_clock::now();
if (_hypothesis.pe_v.empty()) _hypothesis.pe_v.resize(DetectorSpecs::GetME().NOpDets(), 0.);
if (_hypothesis.pe_v.size() != DetectorSpecs::GetME().NOpDets()) {
if (_hypothesis.pe_v.empty()) _hypothesis.pe_v.resize(DetectorSpecs::GetME()->NOpDets(), 0.);
if (_hypothesis.pe_v.size() != DetectorSpecs::GetME()->NOpDets()) {
throw OpT0FinderException("Hypothesis vector length != PMT count");
}

Expand Down Expand Up @@ -526,7 +524,7 @@ FlashMatch_t QLLMatch::OnePMTMatch(const Flash_t& flash) {
double QLLMatch::CallMinuit(const QCluster_t &tpc, const Flash_t &pmt, const bool init_x0) {

if (_measurement.pe_v.empty()) {
_measurement.pe_v.resize(DetectorSpecs::GetME().NOpDets(), 0.);
_measurement.pe_v.resize(DetectorSpecs::GetME()->NOpDets(), 0.);
}
if (_measurement.pe_v.size() != pmt.pe_v.size()) {
std::cout << _measurement.pe_v.size() << " " << pmt.pe_v.size() << std::endl;
Expand Down Expand Up @@ -572,7 +570,7 @@ FlashMatch_t QLLMatch::OnePMTMatch(const Flash_t& flash) {
double reco_x = (_vol_xmax - _vol_xmin)/2;

if (!init_x0) {
reco_x = (_raw_xmax_pt.x - _raw_xmin_pt.x) / 2. - pmt.time * DetectorSpecs::GetME().DriftVelocity() + _raw_xmin_pt.x;
reco_x = (_raw_xmax_pt.x - _raw_xmin_pt.x) / 2. - pmt.time * DetectorSpecs::GetME()->DriftVelocity() + _raw_xmin_pt.x;

if(reco_x < _vol_xmin || reco_x > _vol_xmax)
return kINVALID_DOUBLE;
Expand Down
10 changes: 5 additions & 5 deletions sbncode/OpT0Finder/flashmatch/Algorithms/QWeightPoint.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,7 @@ namespace flashmatch {
{

if(_vis_array.pe_v.empty())
_vis_array.pe_v.resize(DetectorSpecs::GetME().NOpDets());
_vis_array.pe_v.resize(DetectorSpecs::GetME()->NOpDets());

// Prepare the return values (Mostly QWeightPoint)
FlashMatch_t f;
Expand Down Expand Up @@ -62,10 +62,10 @@ namespace flashmatch {
double vis_pe_sum = _vis_array.TotalPE();

double weighted_z = 0;
for(size_t pmt_index=0; pmt_index<DetectorSpecs::GetME().NOpDets(); ++pmt_index) {
for(size_t pmt_index=0; pmt_index<DetectorSpecs::GetME()->NOpDets(); ++pmt_index) {

if(_vis_array.pe_v[pmt_index]<0) continue;
weighted_z += DetectorSpecs::GetME().PMTPosition(pmt_index)[2] * _vis_array.pe_v[pmt_index] / vis_pe_sum;
weighted_z += DetectorSpecs::GetME()->PMTPosition(pmt_index)[2] * _vis_array.pe_v[pmt_index] / vis_pe_sum;

}

Expand All @@ -81,9 +81,9 @@ namespace flashmatch {

f.tpc_point.x = x_offset;

for(size_t pmt_index=0; pmt_index<DetectorSpecs::GetME().NOpDets(); ++pmt_index) {
for(size_t pmt_index=0; pmt_index<DetectorSpecs::GetME()->NOpDets(); ++pmt_index) {
if(_vis_array.pe_v[pmt_index]<0) continue;
f.tpc_point.y += DetectorSpecs::GetME().PMTPosition(pmt_index)[1] * _vis_array.pe_v[pmt_index] / vis_pe_sum;
f.tpc_point.y += DetectorSpecs::GetME()->PMTPosition(pmt_index)[1] * _vis_array.pe_v[pmt_index] / vis_pe_sum;
}

f.tpc_point.z = weighted_z;
Expand Down
12 changes: 6 additions & 6 deletions sbncode/OpT0Finder/flashmatch/Algorithms/TimeCompatMatch.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -43,13 +43,13 @@ namespace flashmatch {

// Earliest flash time => assume clus_x_max is @ detector X-max boundary
#if USING_LARSOFT == 1
double xmax = DetectorSpecs::GetME().ActiveVolume().Max()[0];
double clus_t_min = (clus_x_max - xmax) / DetectorSpecs::GetME().DriftVelocity();
double clus_t_max = clus_x_min / DetectorSpecs::GetME().DriftVelocity();
double xmax = DetectorSpecs::GetME()->ActiveVolume().Max()[0];
double clus_t_min = (clus_x_max - xmax) / DetectorSpecs::GetME()->DriftVelocity();
double clus_t_max = clus_x_min / DetectorSpecs::GetME()->DriftVelocity();
#else
double xmax = DetectorSpecs::GetME().ActiveVolume().Max()[0];
double clus_t_min = (clus_x_max - xmax) / DetectorSpecs::GetME().DriftVelocity();
double clus_t_max = clus_x_min / DetectorSpecs::GetME().DriftVelocity();
double xmax = DetectorSpecs::GetME()->ActiveVolume().Max()[0];
double clus_t_min = (clus_x_max - xmax) / DetectorSpecs::GetME()->DriftVelocity();
double clus_t_max = clus_x_min / DetectorSpecs::GetME()->DriftVelocity();
#endif
FLASH_INFO() << "Cluster xmax: " << xmax << ", clus_t_min: " << clus_t_min
<< ", clus_t_max: " << clus_t_max
Expand Down
4 changes: 2 additions & 2 deletions sbncode/OpT0Finder/flashmatch/Base/FMWKInterface.h
Original file line number Diff line number Diff line change
Expand Up @@ -41,10 +41,10 @@ namespace flashmatch {
DetectorSpecs(std::string filename="specs.cfg");
~DetectorSpecs(){}

inline static DetectorSpecs& GetME(std::string filename="detector_specs.cfg")
inline static DetectorSpecs* GetME(std::string filename="detector_specs.cfg")
{
if(!_me) _me = new DetectorSpecs(filename);
return *_me;
return _me;
}

/// PMT XYZ position filler
Expand Down