LArSoft  v07_13_02
Liquid Argon Software toolkit - http://larsoft.org/
corner::CornerFinderAlg Class Reference

#include "CornerFinderAlg.h"

Public Member Functions

 CornerFinderAlg (fhicl::ParameterSet const &pset)
 
virtual ~CornerFinderAlg ()
 
void reconfigure (fhicl::ParameterSet const &pset)
 
void GrabWires (std::vector< recob::Wire > const &wireVec, geo::Geometry const &)
 
void get_feature_points (std::vector< recob::EndPoint2D > &, geo::Geometry const &)
 
void get_feature_points_LineIntegralScore (std::vector< recob::EndPoint2D > &, geo::Geometry const &)
 
void get_feature_points_fast (std::vector< recob::EndPoint2D > &, geo::Geometry const &)
 
float line_integral (TH2F const &hist, int x1, float y1, int x2, float y2, float threshold)
 
std::vector< float > line_integrals (trkf::BezierTrack &, size_t Steps, float threshold)
 
TH2F const & GetWireDataHist (unsigned int)
 
TH2F const & GetConversionHist (unsigned int)
 
TH2F const & GetDerivativeXHist (unsigned int)
 
TH2F const & GetDerivativeYHist (unsigned int)
 
TH2D const & GetCornerScoreHist (unsigned int)
 
TH2D const & GetMaxSuppressHist (unsigned int)
 

Private Member Functions

void CleanCornerFinderAlg ()
 
void InitializeGeometry (geo::Geometry const &)
 
void create_image_histo (TH2F const &h_wire_data, TH2F &h_conversion)
 
void create_derivative_histograms (TH2F const &h_conversion, TH2F &h_derivative_x, TH2F &h_derivative_y)
 
void create_cornerScore_histogram (TH2F const &h_derivative_x, TH2F const &h_derivative_y, TH2D &h_cornerScore)
 
size_t perform_maximum_suppression (TH2D const &h_cornerScore, std::vector< recob::EndPoint2D > &corner_vector, std::vector< geo::WireID > wireIDs, geo::View_t view, TH2D &h_maxSuppress, int startx=0, int starty=0)
 
size_t calculate_line_integral_score (TH2F const &h_wire_data, std::vector< recob::EndPoint2D > const &corner_vector, std::vector< recob::EndPoint2D > &corner_lineIntegralScore_vector, TH2F &h_lineIntegralScore)
 
void attach_feature_points (TH2F const &h_wire_data, std::vector< geo::WireID > wireIDs, geo::View_t view, std::vector< recob::EndPoint2D > &, int startx=0, int starty=0)
 
void attach_feature_points_LineIntegralScore (TH2F const &h_wire_data, std::vector< geo::WireID > wireIDs, geo::View_t view, std::vector< recob::EndPoint2D > &)
 
void create_smaller_histos (geo::Geometry const &)
 
void remove_duplicates (std::vector< recob::EndPoint2D > &)
 

Private Attributes

std::string fCalDataModuleLabel
 
std::string fConversion_algorithm
 
std::string fConversion_func
 
float fTrimming_threshold
 
int fTrimming_buffer
 
double fTrimming_totalThreshold
 
int fConversion_func_neighborhood
 
float fConversion_threshold
 
int fConversion_bins_per_input_x
 
int fConversion_bins_per_input_y
 
std::string fDerivative_method
 
int fDerivative_neighborhood
 
std::string fDerivative_BlurFunc
 
int fDerivative_BlurNeighborhood
 
int fCornerScore_neighborhood
 
std::string fCornerScore_algorithm
 
float fCornerScore_Noble_epsilon
 
float fCornerScore_Harris_kappa
 
int fMaxSuppress_neighborhood
 
int fMaxSuppress_threshold
 
float fIntegral_bin_threshold
 
float fIntegral_fraction_threshold
 
std::vector< TH2F > WireData_histos
 
std::vector< TH1D > WireData_histos_ProjectionX
 
std::vector< TH1D > WireData_histos_ProjectionY
 
std::vector< std::tuple< int, TH2F, int, int > > WireData_trimmed_histos
 
std::vector< std::vector< geo::WireID > > WireData_IDs
 
std::vector< TH2F > fConversion_histos
 
std::vector< TH2F > fDerivativeX_histos
 
std::vector< TH2F > fDerivativeY_histos
 
std::vector< TH2D > fCornerScore_histos
 
std::vector< TH2D > fMaxSuppress_histos
 
unsigned int event_number
 
unsigned int run_number
 

Detailed Description

Definition at line 29 of file CornerFinderAlg.h.

Constructor & Destructor Documentation

corner::CornerFinderAlg::CornerFinderAlg ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 45 of file CornerFinderAlg.cxx.

References reconfigure().

46 {
47  this->reconfigure(pset);
48 }
void reconfigure(fhicl::ParameterSet const &pset)
corner::CornerFinderAlg::~CornerFinderAlg ( )
virtual

Definition at line 50 of file CornerFinderAlg.cxx.

References CleanCornerFinderAlg().

50  {
51  this->CleanCornerFinderAlg();
52 }

Member Function Documentation

void corner::CornerFinderAlg::attach_feature_points ( TH2F const &  h_wire_data,
std::vector< geo::WireID wireIDs,
geo::View_t  view,
std::vector< recob::EndPoint2D > &  corner_vector,
int  startx = 0,
int  starty = 0 
)
private

Definition at line 509 of file CornerFinderAlg.cxx.

References create_cornerScore_histogram(), create_derivative_histograms(), create_image_histo(), event_number, fConversion_bins_per_input_x, fConversion_bins_per_input_y, perform_maximum_suppression(), run_number, x_max, and x_min.

Referenced by get_feature_points(), and get_feature_points_fast().

514  {
515 
516 
517  const int x_bins = h_wire_data.GetNbinsX();
518  const float x_min = h_wire_data.GetXaxis()->GetBinLowEdge(1);
519  const float x_max = h_wire_data.GetXaxis()->GetBinUpEdge(x_bins);
520 
521  const int y_bins = h_wire_data.GetNbinsY();
522  const float y_min = h_wire_data.GetYaxis()->GetBinLowEdge(1);
523  const float y_max = h_wire_data.GetYaxis()->GetBinUpEdge(y_bins);
524 
525  const int converted_y_bins = y_bins/fConversion_bins_per_input_y;
526  const int converted_x_bins = x_bins/fConversion_bins_per_input_x;
527 
528  std::stringstream conversion_name; conversion_name << "h_conversion_" << view << "_" << run_number << "_" << event_number;
529  std::stringstream dx_name; dx_name << "h_derivative_x_" << view << "_" << run_number << "_" << event_number;
530  std::stringstream dy_name; dy_name << "h_derivative_y_" << view << "_" << run_number << "_" << event_number;
531  std::stringstream cornerScore_name; cornerScore_name << "h_cornerScore_" << view << "_" << run_number << "_" << event_number;
532  std::stringstream maxSuppress_name; maxSuppress_name << "h_maxSuppress_" << view << "_" << run_number << "_" << event_number;
533 
534  TH2F conversion_histo(conversion_name.str().c_str(),"Image Conversion Histogram",
535  converted_x_bins,x_min,x_max,
536  converted_y_bins,y_min,y_max);
537 
538  TH2F derivativeX_histo(dx_name.str().c_str(),"Partial Derivatives (x)",
539  converted_x_bins,x_min,x_max,
540  converted_y_bins,y_min,y_max);
541 
542  TH2F derivativeY_histo(dy_name.str().c_str(),"Partial Derivatives (y)",
543  converted_x_bins,x_min,x_max,
544  converted_y_bins,y_min,y_max);
545 
546  TH2D cornerScore_histo(cornerScore_name.str().c_str(),"Corner Score",
547  converted_x_bins,x_min,x_max,
548  converted_y_bins,y_min,y_max);
549 
550  TH2D maxSuppress_histo(maxSuppress_name.str().c_str(),"Corner Points (Maximum Suppressed)",
551  converted_x_bins,x_min,x_max,
552  converted_y_bins,y_min,y_max);
553 
554  create_image_histo(h_wire_data,conversion_histo);
555  create_derivative_histograms(conversion_histo,derivativeX_histo,derivativeY_histo);
556  create_cornerScore_histogram(derivativeX_histo,derivativeY_histo,cornerScore_histo);
557  perform_maximum_suppression(cornerScore_histo,corner_vector,wireIDs,view,maxSuppress_histo,startx,starty);
558 }
void create_cornerScore_histogram(TH2F const &h_derivative_x, TH2F const &h_derivative_y, TH2D &h_cornerScore)
void create_derivative_histograms(TH2F const &h_conversion, TH2F &h_derivative_x, TH2F &h_derivative_y)
double x_min
Definition: berger.C:15
size_t perform_maximum_suppression(TH2D const &h_cornerScore, std::vector< recob::EndPoint2D > &corner_vector, std::vector< geo::WireID > wireIDs, geo::View_t view, TH2D &h_maxSuppress, int startx=0, int starty=0)
double x_max
Definition: berger.C:16
void create_image_histo(TH2F const &h_wire_data, TH2F &h_conversion)
void corner::CornerFinderAlg::attach_feature_points_LineIntegralScore ( TH2F const &  h_wire_data,
std::vector< geo::WireID wireIDs,
geo::View_t  view,
std::vector< recob::EndPoint2D > &  corner_vector 
)
private

Definition at line 563 of file CornerFinderAlg.cxx.

References calculate_line_integral_score(), create_cornerScore_histogram(), create_derivative_histograms(), create_image_histo(), event_number, fConversion_bins_per_input_x, fConversion_bins_per_input_y, perform_maximum_suppression(), run_number, x_max, and x_min.

Referenced by get_feature_points_LineIntegralScore().

566  {
567 
568 
569  const int x_bins = h_wire_data.GetNbinsX();
570  const float x_min = h_wire_data.GetXaxis()->GetBinLowEdge(1);
571  const float x_max = h_wire_data.GetXaxis()->GetBinUpEdge(x_bins);
572 
573  const int y_bins = h_wire_data.GetNbinsY();
574  const float y_min = h_wire_data.GetYaxis()->GetBinLowEdge(1);
575  const float y_max = h_wire_data.GetYaxis()->GetBinUpEdge(y_bins);
576 
577  const int converted_y_bins = y_bins/fConversion_bins_per_input_y;
578  const int converted_x_bins = x_bins/fConversion_bins_per_input_x;
579 
580  std::stringstream conversion_name; conversion_name << "h_conversion_" << view << "_" << run_number << "_" << event_number;
581  std::stringstream dx_name; dx_name << "h_derivative_x_" << view << "_" << run_number << "_" << event_number;
582  std::stringstream dy_name; dy_name << "h_derivative_y_" << view << "_" << run_number << "_" << event_number;
583  std::stringstream cornerScore_name; cornerScore_name << "h_cornerScore_" << view << "_" << run_number << "_" << event_number;
584  std::stringstream maxSuppress_name; maxSuppress_name << "h_maxSuppress_" << view << "_" << run_number << "_" << event_number;
585 
586  TH2F h_conversion ((conversion_name.str()).c_str(),
587  "Image Conversion Histogram",
588  converted_x_bins,x_min,x_max,
589  converted_y_bins,y_min,y_max);
590  TH2F h_derivative_x((dx_name.str()).c_str(),
591  "Partial Derivatives (x)",
592  converted_x_bins,x_min,x_max,
593  converted_y_bins,y_min,y_max);
594  TH2F h_derivative_y((dy_name.str()).c_str(),
595  "Partial Derivatives (y)",
596  converted_x_bins,x_min,x_max,
597  converted_y_bins,y_min,y_max);
598  TH2D h_cornerScore ((cornerScore_name.str()).c_str(),
599  "Feature Point Corner Score",
600  converted_x_bins,x_min,x_max,
601  converted_y_bins,y_min,y_max);
602  TH2D h_maxSuppress ((maxSuppress_name.str()).c_str(),
603  "Corner Points (Maximum Suppressed)",
604  converted_x_bins,x_min,x_max,
605  converted_y_bins,y_min,y_max);
606 
607  create_image_histo(h_wire_data,h_conversion);
608  create_derivative_histograms(h_conversion,h_derivative_x,h_derivative_y);
609  create_cornerScore_histogram(h_derivative_x,h_derivative_y,h_cornerScore);
610 
611  std::vector<recob::EndPoint2D> corner_vector_tmp;
612  perform_maximum_suppression(h_cornerScore,corner_vector_tmp,wireIDs,view,h_maxSuppress);
613 
614  std::stringstream LI_name; LI_name << "h_lineIntegralScore_" << view << "_" << run_number << "_" << event_number;
615  TH2F h_lineIntegralScore((LI_name.str()).c_str(),
616  "Line Integral Score",
617  x_bins,x_min,x_max,
618  y_bins,y_min,y_max);
619  calculate_line_integral_score(h_wire_data,corner_vector_tmp,corner_vector,h_lineIntegralScore);
620 
621 }
void create_cornerScore_histogram(TH2F const &h_derivative_x, TH2F const &h_derivative_y, TH2D &h_cornerScore)
void create_derivative_histograms(TH2F const &h_conversion, TH2F &h_derivative_x, TH2F &h_derivative_y)
double x_min
Definition: berger.C:15
size_t perform_maximum_suppression(TH2D const &h_cornerScore, std::vector< recob::EndPoint2D > &corner_vector, std::vector< geo::WireID > wireIDs, geo::View_t view, TH2D &h_maxSuppress, int startx=0, int starty=0)
size_t calculate_line_integral_score(TH2F const &h_wire_data, std::vector< recob::EndPoint2D > const &corner_vector, std::vector< recob::EndPoint2D > &corner_lineIntegralScore_vector, TH2F &h_lineIntegralScore)
double x_max
Definition: berger.C:16
void create_image_histo(TH2F const &h_wire_data, TH2F &h_conversion)
size_t corner::CornerFinderAlg::calculate_line_integral_score ( TH2F const &  h_wire_data,
std::vector< recob::EndPoint2D > const &  corner_vector,
std::vector< recob::EndPoint2D > &  corner_lineIntegralScore_vector,
TH2F &  h_lineIntegralScore 
)
private

Definition at line 1152 of file CornerFinderAlg.cxx.

References fIntegral_bin_threshold, fIntegral_fraction_threshold, and line_integral().

Referenced by attach_feature_points_LineIntegralScore().

1155  {
1156 
1157  float score;
1158 
1159  for(auto const i_corner : corner_vector){
1160 
1161  score=0;
1162 
1163  for(auto const j_corner : corner_vector){
1164 
1165 
1166  if( line_integral(h_wire_data,
1167  i_corner.WireID().Wire,i_corner.DriftTime(),
1168  j_corner.WireID().Wire,j_corner.DriftTime(),
1170  {
1171  score+=1.;
1172  }
1173 
1174  }
1175 
1176  recob::EndPoint2D corner(i_corner.DriftTime(),
1177  i_corner.WireID(),
1178  score,
1179  i_corner.ID(),
1180  i_corner.View(),
1181  i_corner.Charge());
1182 
1183  corner_lineIntegralScore_vector.push_back(corner);
1184 
1185 
1186  h_lineIntegralScore.SetBinContent(h_wire_data.GetXaxis()->FindBin(i_corner.WireID().Wire),
1187  h_wire_data.GetYaxis()->FindBin(i_corner.DriftTime()),
1188  score);
1189 
1190  }
1191 
1192  return corner_lineIntegralScore_vector.size();
1193 }
float line_integral(TH2F const &hist, int x1, float y1, int x2, float y2, float threshold)
void corner::CornerFinderAlg::CleanCornerFinderAlg ( )
private

Definition at line 55 of file CornerFinderAlg.cxx.

References WireData_histos, WireData_histos_ProjectionX, WireData_histos_ProjectionY, WireData_IDs, and WireData_trimmed_histos.

Referenced by InitializeGeometry(), and ~CornerFinderAlg().

56 {
57 
58  //for (auto wd_histo : WireData_histos) delete wd_histo;
59  //for (auto wd_histo : WireData_histos_ProjectionX) delete wd_histo;
60  //for (auto wd_histo : WireData_histos_ProjectionY) delete wd_histo;
61  //for (auto histo : fConversion_histos) delete histo;
62  //for (auto histo : fDerivativeX_histos) delete histo;
63  //for (auto histo : fDerivativeY_histos) delete histo;
64  //for (auto histo : fCornerScore_histos) delete histo;
65  //for (auto histo : fMaxSuppress_histos) delete histo;
66 
67  WireData_histos.clear();
70  WireData_IDs.clear();
71 
72  //for (auto wd_histo : WireData_trimmed_histos) delete std::get<1>(wd_histo);
74 
75 }
std::vector< std::tuple< int, TH2F, int, int > > WireData_trimmed_histos
std::vector< TH1D > WireData_histos_ProjectionY
std::vector< TH2F > WireData_histos
std::vector< std::vector< geo::WireID > > WireData_IDs
std::vector< TH1D > WireData_histos_ProjectionX
void corner::CornerFinderAlg::create_cornerScore_histogram ( TH2F const &  h_derivative_x,
TH2F const &  h_derivative_y,
TH2D &  h_cornerScore 
)
private

Definition at line 931 of file CornerFinderAlg.cxx.

References fCornerScore_algorithm, fCornerScore_Harris_kappa, fCornerScore_neighborhood, and fCornerScore_Noble_epsilon.

Referenced by attach_feature_points(), and attach_feature_points_LineIntegralScore().

931  {
932 
933  const int x_bins = h_derivative_x.GetNbinsX();
934  const int y_bins = h_derivative_y.GetNbinsY();
935 
936  //the structure tensor elements
937  double st_xx = 0., st_xy = 0., st_yy = 0.;
938 
939  for(int iy=1+fCornerScore_neighborhood; iy<=(y_bins-fCornerScore_neighborhood); iy++){
940  for(int ix=1+fCornerScore_neighborhood; ix<=(x_bins-fCornerScore_neighborhood); ix++){
941 
942  if(ix==1+fCornerScore_neighborhood){
943  st_xx=0.; st_xy=0.; st_yy=0.;
944 
945  for(int jx=ix-fCornerScore_neighborhood; jx<=ix+fCornerScore_neighborhood; jx++){
946  for(int jy=iy-fCornerScore_neighborhood; jy<=iy+fCornerScore_neighborhood; jy++){
947 
948  st_xx += h_derivative_x.GetBinContent(jx,jy)*h_derivative_x.GetBinContent(jx,jy);
949  st_yy += h_derivative_y.GetBinContent(jx,jy)*h_derivative_y.GetBinContent(jx,jy);
950  st_xy += h_derivative_x.GetBinContent(jx,jy)*h_derivative_y.GetBinContent(jx,jy);
951 
952  }
953  }
954  }
955 
956  // we do it this way to reduce computation time
957  else{
958  for(int jy=iy-fCornerScore_neighborhood; jy<=iy+fCornerScore_neighborhood; jy++){
959 
960  st_xx -= h_derivative_x.GetBinContent(ix-fCornerScore_neighborhood-1,jy)*h_derivative_x.GetBinContent(ix-fCornerScore_neighborhood-1,jy);
961  st_xx += h_derivative_x.GetBinContent(ix+fCornerScore_neighborhood,jy)*h_derivative_x.GetBinContent(ix+fCornerScore_neighborhood,jy);
962 
963  st_yy -= h_derivative_y.GetBinContent(ix-fCornerScore_neighborhood-1,jy)*h_derivative_y.GetBinContent(ix-fCornerScore_neighborhood-1,jy);
964  st_yy += h_derivative_y.GetBinContent(ix+fCornerScore_neighborhood,jy)*h_derivative_y.GetBinContent(ix+fCornerScore_neighborhood,jy);
965 
966  st_xy -= h_derivative_x.GetBinContent(ix-fCornerScore_neighborhood-1,jy)*h_derivative_y.GetBinContent(ix-fCornerScore_neighborhood-1,jy);
967  st_xy += h_derivative_x.GetBinContent(ix+fCornerScore_neighborhood,jy)*h_derivative_y.GetBinContent(ix+fCornerScore_neighborhood,jy);
968  }
969  }
970 
971  if( fCornerScore_algorithm.compare("Noble")==0 ) {
972  h_cornerScore.SetBinContent(ix,iy,
973  (st_xx*st_yy-st_xy*st_xy) / (st_xx+st_yy + fCornerScore_Noble_epsilon));
974  }
975  else if( fCornerScore_algorithm.compare("Harris")==0 ) {
976  h_cornerScore.SetBinContent(ix,iy,
977  (st_xx*st_yy-st_xy*st_xy) - ((st_xx+st_yy)*(st_xx+st_yy)*fCornerScore_Harris_kappa));
978  }
979  else{
980  mf::LogError("CornerFinderAlg") << "BAD CORNER ALGORITHM: " << fCornerScore_algorithm;
981  return;
982  }
983 
984  } // end for loop over x bins
985  } // end for loop over y bins
986 
987 }
std::string fCornerScore_algorithm
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
void corner::CornerFinderAlg::create_derivative_histograms ( TH2F const &  h_conversion,
TH2F &  h_derivative_x,
TH2F &  h_derivative_y 
)
private

Definition at line 686 of file CornerFinderAlg.cxx.

References fDerivative_BlurNeighborhood, fDerivative_method, and fDerivative_neighborhood.

Referenced by attach_feature_points(), and attach_feature_points_LineIntegralScore().

686  {
687 
688  const int x_bins = h_conversion.GetNbinsX();
689  const int y_bins = h_conversion.GetNbinsY();
690 
691  for(int iy=1+fDerivative_neighborhood; iy<=(y_bins-fDerivative_neighborhood); iy++){
692  for(int ix=1+fDerivative_neighborhood; ix<=(x_bins-fDerivative_neighborhood); ix++){
693 
694  if(fDerivative_method.compare("Sobel")==0){
695 
697  h_derivative_x.SetBinContent(ix,iy,
698  0.5*(h_conversion.GetBinContent(ix+1,iy)-h_conversion.GetBinContent(ix-1,iy))
699  + 0.25*(h_conversion.GetBinContent(ix+1,iy+1)-h_conversion.GetBinContent(ix-1,iy+1))
700  + 0.25*(h_conversion.GetBinContent(ix+1,iy-1)-h_conversion.GetBinContent(ix-1,iy-1)));
701  h_derivative_y.SetBinContent(ix,iy,
702  0.5*(h_conversion.GetBinContent(ix,iy+1)-h_conversion.GetBinContent(ix,iy-1))
703  + 0.25*(h_conversion.GetBinContent(ix-1,iy+1)-h_conversion.GetBinContent(ix-1,iy-1))
704  + 0.25*(h_conversion.GetBinContent(ix+1,iy+1)-h_conversion.GetBinContent(ix+1,iy-1)));
705  }
706  else if(fDerivative_neighborhood==2){
707  h_derivative_x.SetBinContent(ix,iy,
708  12*(h_conversion.GetBinContent(ix+1,iy)-h_conversion.GetBinContent(ix-1,iy))
709  + 8*(h_conversion.GetBinContent(ix+1,iy+1)-h_conversion.GetBinContent(ix-1,iy+1))
710  + 8*(h_conversion.GetBinContent(ix+1,iy-1)-h_conversion.GetBinContent(ix-1,iy-1))
711  + 2*(h_conversion.GetBinContent(ix+1,iy+2)-h_conversion.GetBinContent(ix-1,iy+2))
712  + 2*(h_conversion.GetBinContent(ix+1,iy-2)-h_conversion.GetBinContent(ix-1,iy-2))
713  + 6*(h_conversion.GetBinContent(ix+2,iy)-h_conversion.GetBinContent(ix-2,iy))
714  + 4*(h_conversion.GetBinContent(ix+2,iy+1)-h_conversion.GetBinContent(ix-2,iy+1))
715  + 4*(h_conversion.GetBinContent(ix+2,iy-1)-h_conversion.GetBinContent(ix-2,iy-1))
716  + 1*(h_conversion.GetBinContent(ix+2,iy+2)-h_conversion.GetBinContent(ix-2,iy+2))
717  + 1*(h_conversion.GetBinContent(ix+2,iy-2)-h_conversion.GetBinContent(ix-2,iy-2)));
718  h_derivative_y.SetBinContent(ix,iy,
719  12*(h_conversion.GetBinContent(ix,iy+1)-h_conversion.GetBinContent(ix,iy-1))
720  + 8*(h_conversion.GetBinContent(ix-1,iy+1)-h_conversion.GetBinContent(ix-1,iy-1))
721  + 8*(h_conversion.GetBinContent(ix+1,iy+1)-h_conversion.GetBinContent(ix+1,iy-1))
722  + 2*(h_conversion.GetBinContent(ix-2,iy+1)-h_conversion.GetBinContent(ix-2,iy-1))
723  + 2*(h_conversion.GetBinContent(ix+2,iy+1)-h_conversion.GetBinContent(ix+2,iy-1))
724  + 6*(h_conversion.GetBinContent(ix,iy+2)-h_conversion.GetBinContent(ix,iy-2))
725  + 4*(h_conversion.GetBinContent(ix-1,iy+2)-h_conversion.GetBinContent(ix-1,iy-2))
726  + 4*(h_conversion.GetBinContent(ix+1,iy+2)-h_conversion.GetBinContent(ix+1,iy-2))
727  + 1*(h_conversion.GetBinContent(ix-2,iy+2)-h_conversion.GetBinContent(ix-2,iy-2))
728  + 1*(h_conversion.GetBinContent(ix+2,iy+2)-h_conversion.GetBinContent(ix+2,iy-2)));
729  }
730  else{
731  mf::LogError("CornerFinderAlg") << "Sobel derivative not supported for neighborhoods > 2.";
732  return;
733  }
734 
735  } //end if Sobel
736 
737  else if(fDerivative_method.compare("local")==0){
738 
740  h_derivative_x.SetBinContent(ix,iy,
741  (h_conversion.GetBinContent(ix+1,iy)-h_conversion.GetBinContent(ix-1,iy)));
742  h_derivative_y.SetBinContent(ix,iy,
743  (h_conversion.GetBinContent(ix,iy+1)-h_conversion.GetBinContent(ix,iy-1)));
744  }
745  else{
746  mf::LogError("CornerFinderAlg") << "Local derivative not yet supported for neighborhoods > 1.";
747  return;
748  }
749  } //end if local
750 
751  else{
752  mf::LogError("CornerFinderAlg") << "Bad derivative algorithm! " << fDerivative_method;
753  return;
754  }
755 
756  }
757  }
758 
759 
760  //this is just a double Gaussian
761  float func_blur[11][11];
762  func_blur[0][0] = 0.000000;
763  func_blur[0][1] = 0.000000;
764  func_blur[0][2] = 0.000000;
765  func_blur[0][3] = 0.000001;
766  func_blur[0][4] = 0.000002;
767  func_blur[0][5] = 0.000004;
768  func_blur[0][6] = 0.000002;
769  func_blur[0][7] = 0.000001;
770  func_blur[0][8] = 0.000000;
771  func_blur[0][9] = 0.000000;
772  func_blur[0][10] = 0.000000;
773  func_blur[1][0] = 0.000000;
774  func_blur[1][1] = 0.000000;
775  func_blur[1][2] = 0.000004;
776  func_blur[1][3] = 0.000045;
777  func_blur[1][4] = 0.000203;
778  func_blur[1][5] = 0.000335;
779  func_blur[1][6] = 0.000203;
780  func_blur[1][7] = 0.000045;
781  func_blur[1][8] = 0.000004;
782  func_blur[1][9] = 0.000000;
783  func_blur[1][10] = 0.000000;
784  func_blur[2][0] = 0.000000;
785  func_blur[2][1] = 0.000004;
786  func_blur[2][2] = 0.000123;
787  func_blur[2][3] = 0.001503;
788  func_blur[2][4] = 0.006738;
789  func_blur[2][5] = 0.011109;
790  func_blur[2][6] = 0.006738;
791  func_blur[2][7] = 0.001503;
792  func_blur[2][8] = 0.000123;
793  func_blur[2][9] = 0.000004;
794  func_blur[2][10] = 0.000000;
795  func_blur[3][0] = 0.000001;
796  func_blur[3][1] = 0.000045;
797  func_blur[3][2] = 0.001503;
798  func_blur[3][3] = 0.018316;
799  func_blur[3][4] = 0.082085;
800  func_blur[3][5] = 0.135335;
801  func_blur[3][6] = 0.082085;
802  func_blur[3][7] = 0.018316;
803  func_blur[3][8] = 0.001503;
804  func_blur[3][9] = 0.000045;
805  func_blur[3][10] = 0.000001;
806  func_blur[4][0] = 0.000002;
807  func_blur[4][1] = 0.000203;
808  func_blur[4][2] = 0.006738;
809  func_blur[4][3] = 0.082085;
810  func_blur[4][4] = 0.367879;
811  func_blur[4][5] = 0.606531;
812  func_blur[4][6] = 0.367879;
813  func_blur[4][7] = 0.082085;
814  func_blur[4][8] = 0.006738;
815  func_blur[4][9] = 0.000203;
816  func_blur[4][10] = 0.000002;
817  func_blur[5][0] = 0.000004;
818  func_blur[5][1] = 0.000335;
819  func_blur[5][2] = 0.011109;
820  func_blur[5][3] = 0.135335;
821  func_blur[5][4] = 0.606531;
822  func_blur[5][5] = 1.000000;
823  func_blur[5][6] = 0.606531;
824  func_blur[5][7] = 0.135335;
825  func_blur[5][8] = 0.011109;
826  func_blur[5][9] = 0.000335;
827  func_blur[5][10] = 0.000004;
828  func_blur[6][0] = 0.000002;
829  func_blur[6][1] = 0.000203;
830  func_blur[6][2] = 0.006738;
831  func_blur[6][3] = 0.082085;
832  func_blur[6][4] = 0.367879;
833  func_blur[6][5] = 0.606531;
834  func_blur[6][6] = 0.367879;
835  func_blur[6][7] = 0.082085;
836  func_blur[6][8] = 0.006738;
837  func_blur[6][9] = 0.000203;
838  func_blur[6][10] = 0.000002;
839  func_blur[7][0] = 0.000001;
840  func_blur[7][1] = 0.000045;
841  func_blur[7][2] = 0.001503;
842  func_blur[7][3] = 0.018316;
843  func_blur[7][4] = 0.082085;
844  func_blur[7][5] = 0.135335;
845  func_blur[7][6] = 0.082085;
846  func_blur[7][7] = 0.018316;
847  func_blur[7][8] = 0.001503;
848  func_blur[7][9] = 0.000045;
849  func_blur[7][10] = 0.000001;
850  func_blur[8][0] = 0.000000;
851  func_blur[8][1] = 0.000004;
852  func_blur[8][2] = 0.000123;
853  func_blur[8][3] = 0.001503;
854  func_blur[8][4] = 0.006738;
855  func_blur[8][5] = 0.011109;
856  func_blur[8][6] = 0.006738;
857  func_blur[8][7] = 0.001503;
858  func_blur[8][8] = 0.000123;
859  func_blur[8][9] = 0.000004;
860  func_blur[8][10] = 0.000000;
861  func_blur[9][0] = 0.000000;
862  func_blur[9][1] = 0.000000;
863  func_blur[9][2] = 0.000004;
864  func_blur[9][3] = 0.000045;
865  func_blur[9][4] = 0.000203;
866  func_blur[9][5] = 0.000335;
867  func_blur[9][6] = 0.000203;
868  func_blur[9][7] = 0.000045;
869  func_blur[9][8] = 0.000004;
870  func_blur[9][9] = 0.000000;
871  func_blur[9][10] = 0.000000;
872  func_blur[10][0] = 0.000000;
873  func_blur[10][1] = 0.000000;
874  func_blur[10][2] = 0.000000;
875  func_blur[10][3] = 0.000001;
876  func_blur[10][4] = 0.000002;
877  func_blur[10][5] = 0.000004;
878  func_blur[10][6] = 0.000002;
879  func_blur[10][7] = 0.000001;
880  func_blur[10][8] = 0.000000;
881  func_blur[10][9] = 0.000000;
882  func_blur[10][10] = 0.000000;
883 
884  double temp_integral_x = 0;
885  double temp_integral_y = 0;
886 
888 
890  mf::LogWarning("CornerFinderAlg") << "WARNING...BlurNeighborhoods>10 not currently allowed. Shrinking to 10.";
892  }
893 
894  TH2F *h_clone_derivative_x = (TH2F*)h_derivative_x.Clone("h_clone_derivative_x");
895  TH2F *h_clone_derivative_y = (TH2F*)h_derivative_y.Clone("h_clone_derivative_y");
896 
897  temp_integral_x = 0;
898  temp_integral_y = 0;
899 
900  for(int ix=1; ix<=h_derivative_x.GetNbinsX(); ix++){
901  for(int iy=1; iy<=h_derivative_y.GetNbinsY(); iy++){
902 
903  temp_integral_x = 0;
904  temp_integral_y = 0;
905 
906  for(int jx=ix-fDerivative_BlurNeighborhood; jx<=ix+fDerivative_BlurNeighborhood; jx++){
907  for(int jy=iy-fDerivative_BlurNeighborhood; jy<=iy+fDerivative_BlurNeighborhood; jy++){
908  temp_integral_x += h_clone_derivative_x->GetBinContent(jx,jy)*func_blur[(ix-jx)+5][(iy-jy)+5];
909  temp_integral_y += h_clone_derivative_y->GetBinContent(jx,jy)*func_blur[(ix-jx)+5][(iy-jy)+5];
910  }
911  }
912  h_derivative_x.SetBinContent(ix,iy,temp_integral_x);
913  h_derivative_y.SetBinContent(ix,iy,temp_integral_y);
914 
915  }
916  }
917 
918  delete h_clone_derivative_x;
919  delete h_clone_derivative_y;
920 
921 
922  } //end if blur
923 
924 
925 }
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
std::string fDerivative_method
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
void corner::CornerFinderAlg::create_image_histo ( TH2F const &  h_wire_data,
TH2F &  h_conversion 
)
private

Definition at line 626 of file CornerFinderAlg.cxx.

References fConversion_algorithm, fConversion_func, fConversion_func_neighborhood, and fConversion_threshold.

Referenced by attach_feature_points(), and attach_feature_points_LineIntegralScore().

626  {
627 
628  double temp_integral=0;
629 
630  const TF2 fConversion_TF2("fConversion_func",fConversion_func.c_str(),-20,20,-20,20);
631 
632  for(int ix=1; ix<=h_conversion.GetNbinsX(); ix++){
633  for(int iy=1; iy<=h_conversion.GetNbinsY(); iy++){
634 
635  temp_integral = h_wire_data.Integral(ix,ix,iy,iy);
636 
637  if( temp_integral > fConversion_threshold){
638 
639  if(fConversion_algorithm.compare("binary")==0)
640  h_conversion.SetBinContent(ix,iy,10*fConversion_threshold);
641  else if(fConversion_algorithm.compare("standard")==0)
642  h_conversion.SetBinContent(ix,iy,temp_integral);
643 
644  else if(fConversion_algorithm.compare("function")==0){
645 
646  temp_integral = 0;
648  for(int jy=iy-fConversion_func_neighborhood; jy<=iy+fConversion_func_neighborhood; jy++){
649  temp_integral += h_wire_data.GetBinContent(jx,jy)*fConversion_TF2.Eval(ix-jx,iy-jy);
650  }
651  }
652  h_conversion.SetBinContent(ix,iy,temp_integral);
653  }
654 
655  else if(fConversion_algorithm.compare("skeleton")==0){
656 
657  if( (temp_integral > h_wire_data.GetBinContent(ix-1,iy) && temp_integral > h_wire_data.GetBinContent(ix+1,iy))
658  || (temp_integral > h_wire_data.GetBinContent(ix,iy-1) && temp_integral > h_wire_data.GetBinContent(ix,iy+1)))
659  h_conversion.SetBinContent(ix,iy,temp_integral);
660  else
661  h_conversion.SetBinContent(ix,iy,fConversion_threshold);
662  }
663  else if(fConversion_algorithm.compare("sk_bin")==0){
664 
665  if( (temp_integral > h_wire_data.GetBinContent(ix-1,iy) && temp_integral > h_wire_data.GetBinContent(ix+1,iy))
666  || (temp_integral > h_wire_data.GetBinContent(ix,iy-1) && temp_integral > h_wire_data.GetBinContent(ix,iy+1)))
667  h_conversion.SetBinContent(ix,iy,10*fConversion_threshold);
668  else
669  h_conversion.SetBinContent(ix,iy,fConversion_threshold);
670  }
671  else
672  h_conversion.SetBinContent(ix,iy,temp_integral);
673  }
674 
675  else
676  h_conversion.SetBinContent(ix,iy,fConversion_threshold);
677 
678  }
679  }
680 
681 }
std::string fConversion_algorithm
void corner::CornerFinderAlg::create_smaller_histos ( geo::Geometry const &  my_geometry)
private

Definition at line 330 of file CornerFinderAlg.cxx.

References fTrimming_buffer, fTrimming_threshold, fTrimming_totalThreshold, geo::GeometryCore::IteratePlaneIDs(), LOG_DEBUG, WireData_histos, WireData_histos_ProjectionX, WireData_histos_ProjectionY, and WireData_trimmed_histos.

Referenced by get_feature_points_fast().

330  {
331 
332  for(auto const& pid : my_geometry.IteratePlaneIDs() ){
333 
334  LOG_DEBUG("CornerFinderAlg")
335  << "Working plane " << pid.Plane << ".";
336 
337  int x_bins = WireData_histos_ProjectionX.at(pid.Plane).GetNbinsX();
338  int y_bins = WireData_histos_ProjectionY.at(pid.Plane).GetNbinsX();
339 
340  std::vector<int> cut_points_x {0};
341  std::vector<int> cut_points_y {0};
342 
343  for (int ix=1; ix<=x_bins; ix++){
344 
345  float this_value = WireData_histos_ProjectionX.at(pid.Plane).GetBinContent(ix);
346 
347  if(ix<fTrimming_buffer || ix>(x_bins-fTrimming_buffer)) continue;
348 
349  int jx=ix-fTrimming_buffer;
350  while(this_value<fTrimming_threshold){
351  if(jx==ix+fTrimming_buffer) break;
352  this_value = WireData_histos_ProjectionX.at(pid.Plane).GetBinContent(jx);
353  jx++;
354  }
355  if(this_value<fTrimming_threshold){
356  cut_points_x.push_back(ix);
357  }
358 
359  }
360 
361  for (int iy=1; iy<=y_bins; iy++){
362 
363  float this_value = WireData_histos_ProjectionY.at(pid.Plane).GetBinContent(iy);
364 
365  if(iy<fTrimming_buffer || iy>(y_bins-fTrimming_buffer)) continue;
366 
367  int jy=iy-fTrimming_buffer;
368  while(this_value<fTrimming_threshold){
369  if(jy==iy+fTrimming_buffer) break;
370  this_value = WireData_histos_ProjectionY.at(pid.Plane).GetBinContent(jy);
371  jy++;
372  }
373  if(this_value<fTrimming_threshold){
374  cut_points_y.push_back(iy);
375  }
376 
377  }
378 
379  LOG_DEBUG("CornerFinderAlg")
380  << "We have a total of " << cut_points_x.size() << " x cut points."
381  << "\nWe have a total of " << cut_points_y.size() << " y cut points.";
382 
383  std::vector<int> x_low{1};
384  std::vector<int> x_high{x_bins};
385  std::vector<int> y_low{1};
386  std::vector<int> y_high{y_bins};
387  bool x_change = true;
388  bool y_change = true;
389  while(x_change || y_change){
390 
391  x_change = false;
392  y_change = false;
393 
394  size_t current_size = x_low.size();
395 
396  for(size_t il=0; il<current_size; il++){
397 
398  int comp_value = (x_high.at(il) + x_low.at(il)) / 2;
399  std::sort(cut_points_x.begin(),cut_points_x.end(),compare_to_value(comp_value));
400 
401  if(cut_points_x.at(0) <= x_low.at(il) || cut_points_x.at(0) >= x_high.at(il))
402  continue;
403 
404  double integral_low = WireData_histos.at(pid.Plane).Integral(x_low.at(il),cut_points_x.at(0),y_low.at(il),y_high.at(il));
405  double integral_high = WireData_histos.at(pid.Plane).Integral(cut_points_x.at(0),x_high.at(il),y_low.at(il),y_high.at(il));
406  if(integral_low > fTrimming_totalThreshold && integral_high > fTrimming_totalThreshold){
407  x_low.push_back(cut_points_x.at(0));
408  x_high.push_back(x_high.at(il));
409  y_low.push_back(y_low.at(il));
410  y_high.push_back(y_high.at(il));
411 
412  x_high[il] = cut_points_x.at(0);
413  x_change = true;
414  }
415  else if(integral_low > fTrimming_totalThreshold && integral_high < fTrimming_totalThreshold){
416  x_high[il] = cut_points_x.at(0);
417  x_change = true;
418  }
419  else if(integral_low < fTrimming_totalThreshold && integral_high > fTrimming_totalThreshold){
420  x_low[il] = cut_points_x.at(0);
421  x_change = true;
422  }
423  }
424 
425  current_size = x_low.size();
426 
427  for(size_t il=0; il<current_size; il++){
428 
429  int comp_value = (y_high.at(il) - y_low.at(il)) / 2;
430  std::sort(cut_points_y.begin(),cut_points_y.end(),compare_to_value(comp_value));
431 
432  if(cut_points_y.at(0) <= y_low.at(il) || cut_points_y.at(0) >= y_high.at(il))
433  continue;
434 
435  double integral_low = WireData_histos.at(pid.Plane).Integral(x_low.at(il),x_high.at(il),y_low.at(il),cut_points_y.at(0));
436  double integral_high = WireData_histos.at(pid.Plane).Integral(x_low.at(il),x_high.at(il),cut_points_y.at(0),y_high.at(il));
437  if(integral_low > fTrimming_totalThreshold && integral_high > fTrimming_totalThreshold){
438  y_low.push_back(cut_points_y.at(0));
439  y_high.push_back(y_high.at(il));
440  x_low.push_back(x_low.at(il));
441  x_high.push_back(x_high.at(il));
442 
443  y_high[il] = cut_points_y.at(0);
444  y_change = true;
445  }
446  else if(integral_low > fTrimming_totalThreshold && integral_high < fTrimming_totalThreshold){
447  y_high[il] = cut_points_y.at(0);
448  y_change = true;
449  }
450  else if(integral_low < fTrimming_totalThreshold && integral_high > fTrimming_totalThreshold){
451  y_low[il] = cut_points_y.at(0);
452  y_change = true;
453  }
454  }
455 
456  }
457 
458  LOG_DEBUG("CornerFinderAlg")
459  << "First point in x is " << cut_points_x.at(0);
460 
461  std::sort(cut_points_x.begin(),cut_points_x.end(),compare_to_value(x_bins/2));
462 
463  LOG_DEBUG("CornerFinderAlg")
464  << "Now the first point in x is " << cut_points_x.at(0);
465 
466  LOG_DEBUG("CornerFinderAlg")
467  << "First point in y is " << cut_points_y.at(0);
468 
469  std::sort(cut_points_y.begin(),cut_points_y.end(),compare_to_value(y_bins/2));
470 
471  LOG_DEBUG("CornerFinderAlg")
472  << "Now the first point in y is " << cut_points_y.at(0);
473 
474  LOG_DEBUG("CornerFinderAlg")
475  << "\nIntegral on the SW side is "
476  << WireData_histos.at(pid.Plane).Integral(1,cut_points_x.at(0),1,cut_points_y.at(0))
477  << "\nIntegral on the SE side is "
478  << WireData_histos.at(pid.Plane).Integral(cut_points_x.at(0),x_bins,1,cut_points_y.at(0))
479  << "\nIntegral on the NW side is "
480  << WireData_histos.at(pid.Plane).Integral(1,cut_points_x.at(0),cut_points_y.at(0),y_bins)
481  << "\nIntegral on the NE side is "
482  << WireData_histos.at(pid.Plane).Integral(cut_points_x.at(0),x_bins,cut_points_y.at(0),y_bins);
483 
484 
485  for(size_t il=0; il<x_low.size(); il++){
486 
487  std::stringstream h_name;
488  h_name << "h_" << pid.Cryostat << "_" << pid.TPC << "_" << pid.Plane << "_sub" << il;
489  TH2F h_tmp( (h_name.str()).c_str(),"",
490  x_high.at(il)-x_low.at(il)+1,x_low.at(il),x_high.at(il),
491  y_high.at(il)-y_low.at(il)+1,y_low.at(il),y_high.at(il));
492 
493  for(int ix=1; ix<=(x_high.at(il)-x_low.at(il)+1); ix++){
494  for(int iy=1; iy<=(y_high.at(il)-y_low.at(il)+1); iy++){
495  h_tmp.SetBinContent(ix,iy,WireData_histos.at(pid.Plane).GetBinContent(x_low.at(il)+(ix-1),y_low.at(il)+(iy-1)));
496  }
497  }
498 
499  WireData_trimmed_histos.push_back(std::make_tuple(pid.Plane,h_tmp,x_low.at(il)-1,y_low.at(il)-1));
500  }
501 
502  }// end loop over PlaneIDs
503 
504 
505 }
std::vector< std::tuple< int, TH2F, int, int > > WireData_trimmed_histos
std::vector< TH1D > WireData_histos_ProjectionY
std::vector< TH2F > WireData_histos
std::vector< TH1D > WireData_histos_ProjectionX
#define LOG_DEBUG(id)
void corner::CornerFinderAlg::get_feature_points ( std::vector< recob::EndPoint2D > &  corner_vector,
geo::Geometry const &  my_geometry 
)

Definition at line 201 of file CornerFinderAlg.cxx.

References attach_feature_points(), geo::GeometryCore::IteratePlaneIDs(), geo::GeometryCore::View(), WireData_histos, and WireData_IDs.

Referenced by trkf::FeatureTracker::produce().

202  {
203 
204 
205 
206  for(auto const& pid : my_geometry.IteratePlaneIDs()){
208  WireData_IDs.at(pid.Plane),
209  my_geometry.View(pid),
210  corner_vector);
211  }
212 
213 }
void attach_feature_points(TH2F const &h_wire_data, std::vector< geo::WireID > wireIDs, geo::View_t view, std::vector< recob::EndPoint2D > &, int startx=0, int starty=0)
std::vector< TH2F > WireData_histos
std::vector< std::vector< geo::WireID > > WireData_IDs
void corner::CornerFinderAlg::get_feature_points_fast ( std::vector< recob::EndPoint2D > &  corner_vector,
geo::Geometry const &  my_geometry 
)

Definition at line 217 of file CornerFinderAlg.cxx.

References attach_feature_points(), create_smaller_histos(), geo::GeometryCore::Cryostat(), LOG_DEBUG, geo::GeometryCore::Ncryostats(), geo::CryostatGeo::NTPC(), geo::TPCGeo::Plane(), geo::CryostatGeo::TPC(), geo::PlaneGeo::View(), WireData_IDs, and WireData_trimmed_histos.

Referenced by vertex::CornerFinder::produce().

218  {
219 
220 
221 
222  create_smaller_histos(my_geometry);
223 
224  for(unsigned int cstat = 0; cstat < my_geometry.Ncryostats(); ++cstat){
225  for(unsigned int tpc = 0; tpc < my_geometry.Cryostat(cstat).NTPC(); ++tpc){
226  for(size_t histos=0; histos!= WireData_trimmed_histos.size(); histos++){
227 
228  int plane = std::get<0>(WireData_trimmed_histos.at(histos));
229  int startx = std::get<2>(WireData_trimmed_histos.at(histos));
230  int starty = std::get<3>(WireData_trimmed_histos.at(histos));
231 
232  LOG_DEBUG("CornerFinderAlg")
233  << "Doing histogram " << histos
234  << ", of plane " << plane
235  << " with start points " << startx << " " << starty;
236 
237  attach_feature_points(std::get<1>(WireData_trimmed_histos.at(histos)),
238  WireData_IDs.at(plane),my_geometry.Cryostat(cstat).TPC(tpc).Plane(plane).View(),corner_vector,startx,starty);
239 
240  LOG_DEBUG("CornerFinderAlg") << "Total feature points now is " << corner_vector.size();
241  }
242 
243  //remove_duplicates(corner_vector);
244 
245  }
246  }
247 
248 }
std::vector< std::tuple< int, TH2F, int, int > > WireData_trimmed_histos
void create_smaller_histos(geo::Geometry const &)
void attach_feature_points(TH2F const &h_wire_data, std::vector< geo::WireID > wireIDs, geo::View_t view, std::vector< recob::EndPoint2D > &, int startx=0, int starty=0)
std::vector< std::vector< geo::WireID > > WireData_IDs
#define LOG_DEBUG(id)
void corner::CornerFinderAlg::get_feature_points_LineIntegralScore ( std::vector< recob::EndPoint2D > &  corner_vector,
geo::Geometry const &  my_geometry 
)

Definition at line 253 of file CornerFinderAlg.cxx.

References attach_feature_points_LineIntegralScore(), geo::GeometryCore::IteratePlaneIDs(), geo::GeometryCore::View(), WireData_histos, and WireData_IDs.

254  {
255 
256 
257 
258  for(auto const& pid : my_geometry.IteratePlaneIDs()){
260  WireData_IDs.at(pid.Plane),
261  my_geometry.View(pid),
262  corner_vector);
263  }
264 
265 }
void attach_feature_points_LineIntegralScore(TH2F const &h_wire_data, std::vector< geo::WireID > wireIDs, geo::View_t view, std::vector< recob::EndPoint2D > &)
std::vector< TH2F > WireData_histos
std::vector< std::vector< geo::WireID > > WireData_IDs
TH2F const& corner::CornerFinderAlg::GetConversionHist ( unsigned  int)
TH2D const& corner::CornerFinderAlg::GetCornerScoreHist ( unsigned  int)
TH2F const& corner::CornerFinderAlg::GetDerivativeXHist ( unsigned  int)
TH2F const& corner::CornerFinderAlg::GetDerivativeYHist ( unsigned  int)
TH2D const& corner::CornerFinderAlg::GetMaxSuppressHist ( unsigned  int)
TH2F const & corner::CornerFinderAlg::GetWireDataHist ( unsigned int  i_plane)

Definition at line 1197 of file CornerFinderAlg.cxx.

References WireData_histos.

Referenced by trkf::FeatureTracker::CheckSeedLineInt(), and trkf::FeatureTracker::GetValidLines().

1197  {
1198  return WireData_histos.at(i_plane);
1199 }
std::vector< TH2F > WireData_histos
void corner::CornerFinderAlg::GrabWires ( std::vector< recob::Wire > const &  wireVec,
geo::Geometry const &  my_geometry 
)

Definition at line 139 of file CornerFinderAlg.cxx.

References geo::GeometryCore::ChannelToWire(), fCalDataModuleLabel, InitializeGeometry(), geo::GeometryCore::Nplanes(), geo::GeometryCore::Nwires(), geo::PlaneID::Plane, geo::WireID::Wire, WireData_histos, WireData_histos_ProjectionX, WireData_histos_ProjectionY, and WireData_IDs.

Referenced by vertex::CornerFinder::produce(), and trkf::FeatureTracker::produce().

139  {
140 
141  InitializeGeometry(my_geometry);
142 
143  const unsigned int nTimeTicks = wireVec.at(0).NSignal();
144 
145  // Initialize the histograms.
146  // All of this should eventually be changed to not need to use histograms...
147  for (unsigned int i_plane=0; i_plane < my_geometry.Nplanes(); i_plane++){
148 
149  std::stringstream ss_tmp_name,ss_tmp_title;
150  ss_tmp_name << "h_WireData_" << i_plane;
151  ss_tmp_title << fCalDataModuleLabel << " wire data for plane " << i_plane << ";Wire Number;Time Tick";
152 
153  if( (unsigned int)(WireData_histos.at(i_plane).GetNbinsX()) == (my_geometry.Nwires(i_plane)) ) {
154  WireData_histos.at(i_plane).Reset();
155  WireData_histos.at(i_plane).SetName(ss_tmp_name.str().c_str());
156  WireData_histos.at(i_plane).SetTitle(ss_tmp_title.str().c_str());
157  }
158  else
159  WireData_histos.at(i_plane) = TH2F(ss_tmp_name.str().c_str(),
160  ss_tmp_title.str().c_str(),
161  my_geometry.Nwires(i_plane),
162  0,
163  my_geometry.Nwires(i_plane),
164  nTimeTicks,
165  0,
166  nTimeTicks);
167  }
168 
169 
170  /* Now do the loop over the wires. */
171  for (std::vector<recob::Wire>::const_iterator iwire = wireVec.begin(); iwire < wireVec.end(); iwire++) {
172 
173  std::vector<geo::WireID> possible_wireIDs = my_geometry.ChannelToWire(iwire->Channel());
174  geo::WireID this_wireID;
175  try { this_wireID = possible_wireIDs.at(0);}
176  catch(cet::exception& excep) { mf::LogError("CornerFinderAlg") << "Bail out! No Possible Wires!\n"; }
177 
178  unsigned int i_plane = this_wireID.Plane;
179  unsigned int i_wire = this_wireID.Wire;
180 
181  WireData_IDs.at(i_plane).at(i_wire) = this_wireID;
182 
183  for(unsigned int i_time = 0; i_time < nTimeTicks; i_time++){
184  WireData_histos.at(i_plane).SetBinContent(i_wire,i_time,(iwire->Signal()).at(i_time));
185  }//<---End time loop
186 
187  }//<-- End loop over wires
188 
189 
190  for (unsigned int i_plane=0; i_plane < my_geometry.Nplanes(); i_plane++){
191  WireData_histos_ProjectionX.at(i_plane) = *(WireData_histos.at(i_plane).ProjectionX());
192  WireData_histos_ProjectionY.at(i_plane) = *(WireData_histos.at(i_plane).ProjectionY());
193  }
194 
195 
196 }
std::string fCalDataModuleLabel
void InitializeGeometry(geo::Geometry const &)
WireID_t Wire
Index of the wire within its plane.
Definition: geo_types.h:313
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
intermediate_table::const_iterator const_iterator
std::vector< TH1D > WireData_histos_ProjectionY
std::vector< TH2F > WireData_histos
PlaneID_t Plane
Index of the plane within its TPC.
Definition: geo_types.h:258
std::vector< std::vector< geo::WireID > > WireData_IDs
std::vector< TH1D > WireData_histos_ProjectionX
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void corner::CornerFinderAlg::InitializeGeometry ( geo::Geometry const &  my_geometry)
private

Definition at line 113 of file CornerFinderAlg.cxx.

References CleanCornerFinderAlg(), fConversion_histos, fCornerScore_histos, fDerivativeX_histos, fDerivativeY_histos, fMaxSuppress_histos, geo::GeometryCore::Nplanes(), geo::GeometryCore::Nwires(), WireData_histos, WireData_histos_ProjectionX, WireData_histos_ProjectionY, WireData_IDs, and WireData_trimmed_histos.

Referenced by GrabWires().

113  {
114 
116 
117  // set the sizes of the WireData_histos and WireData_IDs
118  unsigned int nPlanes = my_geometry.Nplanes();
119  WireData_histos.resize(nPlanes);
120  WireData_histos_ProjectionX.resize(nPlanes);
121  WireData_histos_ProjectionY.resize(nPlanes);
122  fConversion_histos.resize(nPlanes);
123  fDerivativeX_histos.resize(nPlanes);
124  fDerivativeY_histos.resize(nPlanes);
125  fCornerScore_histos.resize(nPlanes);
126  fMaxSuppress_histos.resize(nPlanes);
127 
128  /* For now, we need something to associate each wire in the histogram with a wire_id.
129  This is not a beautiful way of handling this, but for now it should work. */
130  WireData_IDs.resize(nPlanes);
131  for(unsigned int i_plane=0; i_plane < nPlanes; ++i_plane)
132  WireData_IDs.at(i_plane).resize(my_geometry.Nwires(i_plane));
133 
134  WireData_trimmed_histos.resize(0);
135 
136 }
std::vector< std::tuple< int, TH2F, int, int > > WireData_trimmed_histos
std::vector< TH2D > fMaxSuppress_histos
std::vector< TH2D > fCornerScore_histos
std::vector< TH2F > fDerivativeX_histos
std::vector< TH2F > fDerivativeY_histos
std::vector< TH1D > WireData_histos_ProjectionY
std::vector< TH2F > WireData_histos
std::vector< std::vector< geo::WireID > > WireData_IDs
std::vector< TH1D > WireData_histos_ProjectionX
std::vector< TH2F > fConversion_histos
float corner::CornerFinderAlg::line_integral ( TH2F const &  hist,
int  x1,
float  y1,
int  x2,
float  y2,
float  threshold 
)

Definition at line 1053 of file CornerFinderAlg.cxx.

References e, tmp, x1, x2, y1, and y2.

Referenced by calculate_line_integral_score(), trkf::FeatureTracker::CheckSeedLineInt(), and trkf::FeatureTracker::GetValidLines().

1053  {
1054 
1055  int x1 = hist.GetXaxis()->FindBin( begin_x );
1056  int y1 = hist.GetYaxis()->FindBin( begin_y );
1057  int x2 = hist.GetXaxis()->FindBin( end_x );
1058  int y2 = hist.GetYaxis()->FindBin( end_y );
1059 
1060  if(x1==x2 && abs(y1-y2)<1e-5)
1061  return 0;
1062 
1063  if(x2<x1){
1064  int tmp = x2;
1065  x2 = x1;
1066  x1 = tmp;
1067 
1068  int tmp_y = y2;
1069  y2 = y1;
1070  y1 = tmp_y;
1071  }
1072 
1073  float fraction = 0;
1074  int bin_counter = 0;
1075 
1076  if(x2!=x1){
1077 
1078  float slope = (y2-y1)/((float)(x2-x1));
1079 
1080  for(int ix=x1; ix<=x2; ix++){
1081 
1082  int y_min,y_max;
1083 
1084  if(slope>=0){
1085  y_min = y1 + slope*(ix-x1);
1086  y_max = y1 + slope*(ix+1-x1);
1087  }
1088  else {
1089  y_max = (y1+1) + slope*(ix-x1);
1090  y_min = (y1+1) + slope*(ix+1-x1);
1091  }
1092 
1093  for(int iy=y_min; iy<=y_max; iy++){
1094  bin_counter++;
1095 
1096  if( hist.GetBinContent(ix,iy) > threshold )
1097  fraction += 1.;
1098  }
1099 
1100  }
1101  }
1102  else{
1103 
1104  int y_min,y_max;
1105  if(y1<y2){
1106  y_min=y1; y_max=y2;
1107  }
1108  else{
1109  y_min=y2; y_max=y1;
1110  }
1111  for(int iy=y_min; iy<=y_max; iy++){
1112  bin_counter++;
1113  if( hist.GetBinContent(x1,iy) > threshold)
1114  fraction += 1.;
1115  }
1116 
1117  }
1118 
1119  return fraction/bin_counter;
1120 }
Float_t y1[n_points_granero]
Definition: compare.C:5
Float_t x1[n_points_granero]
Definition: compare.C:5
Float_t tmp
Definition: plot.C:37
Float_t y2[n_points_geant4]
Definition: compare.C:26
TH2F * hist
Definition: plot.C:136
Float_t x2[n_points_geant4]
Definition: compare.C:26
Float_t e
Definition: plot.C:34
std::vector< float > corner::CornerFinderAlg::line_integrals ( trkf::BezierTrack TheTrack,
size_t  Steps,
float  threshold 
)

Definition at line 1123 of file CornerFinderAlg.cxx.

References trkf::BezierTrack::GetProjectedPointUVWT(), s, WireData_histos, x, and y.

Referenced by trkf::FeatureTracker::GetConnectionMap().

1123  {
1124 
1125  std::vector<float> fractions(3);
1126 
1127  for(size_t i=0; i!=Steps; ++i)
1128  {
1129  double s = float(i)/float(Steps);
1130  double uvw[3], ticks[3];
1131  int c=0, t=0;
1132  TheTrack.GetProjectedPointUVWT(s, uvw, ticks, c, t);
1133 
1134  for(size_t j=0; j!=WireData_histos.size(); ++j)
1135  {
1136  int x = WireData_histos.at(j).GetXaxis()->FindBin(uvw[j]);
1137  int y = WireData_histos.at(j).GetYaxis()->FindBin(ticks[j]);
1138 
1139  if( WireData_histos.at(j).GetBinContent(x,y) > threshold )
1140  fractions.at(j) += 1.;
1141  }
1142  }
1143  for(size_t j=0; j!=fractions.size(); ++j)
1144  fractions.at(j) /= float(Steps);
1145 
1146  return fractions;
1147 }
Float_t x
Definition: compare.C:6
Float_t s
Definition: plot.C:23
Float_t y
Definition: compare.C:6
std::vector< TH2F > WireData_histos
void GetProjectedPointUVWT(double s, double *uvw, double *ticks, int c, int t) const
size_t corner::CornerFinderAlg::perform_maximum_suppression ( TH2D const &  h_cornerScore,
std::vector< recob::EndPoint2D > &  corner_vector,
std::vector< geo::WireID wireIDs,
geo::View_t  view,
TH2D &  h_maxSuppress,
int  startx = 0,
int  starty = 0 
)
private

Definition at line 992 of file CornerFinderAlg.cxx.

References fConversion_bins_per_input_x, fConversion_bins_per_input_y, fMaxSuppress_neighborhood, and fMaxSuppress_threshold.

Referenced by attach_feature_points(), and attach_feature_points_LineIntegralScore().

998  {
999 
1000  const int x_bins = h_cornerScore.GetNbinsX();
1001  const int y_bins = h_cornerScore.GetNbinsY();
1002 
1003  double temp_max;
1004  bool temp_center_bin;
1005 
1006  for(int iy=1; iy<=y_bins; iy++){
1007  for(int ix=1; ix<=x_bins; ix++){
1008 
1009  if(h_cornerScore.GetBinContent(ix,iy) < fMaxSuppress_threshold)
1010  continue;
1011 
1012  temp_max = -1000;
1013  temp_center_bin = false;
1014 
1015  for(int jx=ix-fMaxSuppress_neighborhood; jx<=ix+fMaxSuppress_neighborhood; jx++){
1016  for(int jy=iy-fMaxSuppress_neighborhood; jy<=iy+fMaxSuppress_neighborhood; jy++){
1017 
1018  if(h_cornerScore.GetBinContent(jx,jy) > temp_max){
1019  temp_max = h_cornerScore.GetBinContent(jx,jy);
1020  if(jx==ix && jy==iy) temp_center_bin=true;
1021  else{ temp_center_bin=false; }
1022  }
1023 
1024  }
1025  }
1026 
1027  if(temp_center_bin){
1028 
1029  float time_tick = 0.5 * (float)((2*(iy+starty)) * fConversion_bins_per_input_y);
1030  int wire_number = ( (2*(ix+startx))*fConversion_bins_per_input_x ) / 2;
1031  double totalQ = 0;
1032  int id = 0;
1033  recob::EndPoint2D corner(time_tick,
1034  wireIDs[wire_number],
1035  h_cornerScore.GetBinContent(ix,iy),
1036  id,
1037  view,
1038  totalQ);
1039  corner_vector.push_back(corner);
1040 
1041  h_maxSuppress.SetBinContent(ix,iy,h_cornerScore.GetBinContent(ix,iy));
1042  }
1043 
1044  }
1045  }
1046 
1047  return corner_vector.size();
1048 
1049 }
void corner::CornerFinderAlg::reconfigure ( fhicl::ParameterSet const &  pset)

Definition at line 78 of file CornerFinderAlg.cxx.

References fCalDataModuleLabel, fConversion_algorithm, fConversion_bins_per_input_x, fConversion_bins_per_input_y, fConversion_func, fConversion_func_neighborhood, fConversion_threshold, fCornerScore_algorithm, fCornerScore_Harris_kappa, fCornerScore_neighborhood, fCornerScore_Noble_epsilon, fDerivative_BlurFunc, fDerivative_BlurNeighborhood, fDerivative_method, fDerivative_neighborhood, fIntegral_bin_threshold, fIntegral_fraction_threshold, fMaxSuppress_neighborhood, fMaxSuppress_threshold, fTrimming_buffer, fTrimming_threshold, fTrimming_totalThreshold, and fhicl::ParameterSet::get().

Referenced by CornerFinderAlg(), and vertex::CornerFinder::reconfigure().

79 {
80  // ### These are all the tuneable .fcl file parameters from the event ###
81  fCalDataModuleLabel = p.get< std::string >("CalDataModuleLabel");
82  fTrimming_threshold = p.get< float >("Trimming_threshold");
83  fTrimming_totalThreshold = p.get< double >("Trimming_totalThreshold");
84  fConversion_threshold = p.get< float >("Conversion_threshold");
85  fConversion_bins_per_input_x = p.get< int >("Conversion_bins_per_input_x");
86  fConversion_bins_per_input_y = p.get< int >("Conversion_bins_per_input_y");
87  fConversion_algorithm = p.get< std::string >("Conversion_algorithm");
88  fConversion_func = p.get< std::string >("Conversion_function");
89  fConversion_func_neighborhood = p.get< int >("Conversion_func_neighborhood");
90  fDerivative_method = p.get< std::string >("Derivative_method");
91  fDerivative_neighborhood = p.get< int >("Derivative_neighborhood");
92  fDerivative_BlurFunc = p.get< std::string >("Derivative_BlurFunc");
93  fDerivative_BlurNeighborhood = p.get< int >("Derivative_BlurNeighborhood");
94  fCornerScore_neighborhood = p.get< int >("CornerScore_neighborhood");
95  fCornerScore_algorithm = p.get< std::string >("CornerScore_algorithm");
96  fCornerScore_Noble_epsilon = p.get< float >("CornerScore_Noble_epsilon");
97  fCornerScore_Harris_kappa = p.get< float >("CornerScore_Harris_kappa");
98  fMaxSuppress_neighborhood = p.get< int >("MaxSuppress_neighborhood");
99  fMaxSuppress_threshold = p.get< int >("MaxSuppress_threshold");
100  fIntegral_bin_threshold = p.get< float >("Integral_bin_threshold");
101  fIntegral_fraction_threshold = p.get< float >("Integral_fraction_threshold");
102 
103  int neighborhoods[] = { fConversion_func_neighborhood,
107  fMaxSuppress_neighborhood };
108  fTrimming_buffer = *std::max_element(neighborhoods,neighborhoods+5);
109 
110 }
std::string fDerivative_BlurFunc
std::string fCornerScore_algorithm
std::string fCalDataModuleLabel
std::string fConversion_algorithm
std::string fDerivative_method
void corner::CornerFinderAlg::remove_duplicates ( std::vector< recob::EndPoint2D > &  corner_vector)
private

Definition at line 268 of file CornerFinderAlg.cxx.

268  {
269 
270  int i_wire, j_wire;
271  float i_time, j_time;
272  for(size_t i=0; i != corner_vector.size(); i++){
273 
274  i_wire = corner_vector.at(i).WireID().Wire;
275  i_time = corner_vector.at(i).DriftTime();
276 
277  for(size_t j=i+1; j != corner_vector.size(); j++){
278 
279  j_wire = corner_vector.at(j).WireID().Wire;
280  j_time = corner_vector.at(j).DriftTime();
281 
282  if(std::abs(i_wire-j_wire) < 5 && std::abs(i_time - j_time) < 10){
283  corner_vector.erase(corner_vector.begin()+j);
284  j--;
285  }
286 
287  }
288 
289  }
290 
291 }

Member Data Documentation

unsigned int corner::CornerFinderAlg::event_number
private
std::string corner::CornerFinderAlg::fCalDataModuleLabel
private

Definition at line 70 of file CornerFinderAlg.h.

Referenced by GrabWires(), and reconfigure().

std::string corner::CornerFinderAlg::fConversion_algorithm
private

Definition at line 71 of file CornerFinderAlg.h.

Referenced by create_image_histo(), and reconfigure().

int corner::CornerFinderAlg::fConversion_bins_per_input_x
private
int corner::CornerFinderAlg::fConversion_bins_per_input_y
private
std::string corner::CornerFinderAlg::fConversion_func
private

Definition at line 72 of file CornerFinderAlg.h.

Referenced by create_image_histo(), and reconfigure().

int corner::CornerFinderAlg::fConversion_func_neighborhood
private

Definition at line 76 of file CornerFinderAlg.h.

Referenced by create_image_histo(), and reconfigure().

std::vector<TH2F> corner::CornerFinderAlg::fConversion_histos
private

Definition at line 99 of file CornerFinderAlg.h.

Referenced by InitializeGeometry().

float corner::CornerFinderAlg::fConversion_threshold
private

Definition at line 77 of file CornerFinderAlg.h.

Referenced by create_image_histo(), and reconfigure().

std::string corner::CornerFinderAlg::fCornerScore_algorithm
private

Definition at line 85 of file CornerFinderAlg.h.

Referenced by create_cornerScore_histogram(), and reconfigure().

float corner::CornerFinderAlg::fCornerScore_Harris_kappa
private

Definition at line 87 of file CornerFinderAlg.h.

Referenced by create_cornerScore_histogram(), and reconfigure().

std::vector<TH2D> corner::CornerFinderAlg::fCornerScore_histos
private

Definition at line 102 of file CornerFinderAlg.h.

Referenced by InitializeGeometry().

int corner::CornerFinderAlg::fCornerScore_neighborhood
private

Definition at line 84 of file CornerFinderAlg.h.

Referenced by create_cornerScore_histogram(), and reconfigure().

float corner::CornerFinderAlg::fCornerScore_Noble_epsilon
private

Definition at line 86 of file CornerFinderAlg.h.

Referenced by create_cornerScore_histogram(), and reconfigure().

std::string corner::CornerFinderAlg::fDerivative_BlurFunc
private

Definition at line 82 of file CornerFinderAlg.h.

Referenced by reconfigure().

int corner::CornerFinderAlg::fDerivative_BlurNeighborhood
private

Definition at line 83 of file CornerFinderAlg.h.

Referenced by create_derivative_histograms(), and reconfigure().

std::string corner::CornerFinderAlg::fDerivative_method
private

Definition at line 80 of file CornerFinderAlg.h.

Referenced by create_derivative_histograms(), and reconfigure().

int corner::CornerFinderAlg::fDerivative_neighborhood
private

Definition at line 81 of file CornerFinderAlg.h.

Referenced by create_derivative_histograms(), and reconfigure().

std::vector<TH2F> corner::CornerFinderAlg::fDerivativeX_histos
private

Definition at line 100 of file CornerFinderAlg.h.

Referenced by InitializeGeometry().

std::vector<TH2F> corner::CornerFinderAlg::fDerivativeY_histos
private

Definition at line 101 of file CornerFinderAlg.h.

Referenced by InitializeGeometry().

float corner::CornerFinderAlg::fIntegral_bin_threshold
private

Definition at line 90 of file CornerFinderAlg.h.

Referenced by calculate_line_integral_score(), and reconfigure().

float corner::CornerFinderAlg::fIntegral_fraction_threshold
private

Definition at line 91 of file CornerFinderAlg.h.

Referenced by calculate_line_integral_score(), and reconfigure().

std::vector<TH2D> corner::CornerFinderAlg::fMaxSuppress_histos
private

Definition at line 103 of file CornerFinderAlg.h.

Referenced by InitializeGeometry().

int corner::CornerFinderAlg::fMaxSuppress_neighborhood
private

Definition at line 88 of file CornerFinderAlg.h.

Referenced by perform_maximum_suppression(), and reconfigure().

int corner::CornerFinderAlg::fMaxSuppress_threshold
private

Definition at line 89 of file CornerFinderAlg.h.

Referenced by perform_maximum_suppression(), and reconfigure().

int corner::CornerFinderAlg::fTrimming_buffer
private

Definition at line 74 of file CornerFinderAlg.h.

Referenced by create_smaller_histos(), and reconfigure().

float corner::CornerFinderAlg::fTrimming_threshold
private

Definition at line 73 of file CornerFinderAlg.h.

Referenced by create_smaller_histos(), and reconfigure().

double corner::CornerFinderAlg::fTrimming_totalThreshold
private

Definition at line 75 of file CornerFinderAlg.h.

Referenced by create_smaller_histos(), and reconfigure().

unsigned int corner::CornerFinderAlg::run_number
private
std::vector<TH2F> corner::CornerFinderAlg::WireData_histos
private
std::vector<TH1D> corner::CornerFinderAlg::WireData_histos_ProjectionX
private
std::vector<TH1D> corner::CornerFinderAlg::WireData_histos_ProjectionY
private
std::vector< std::vector<geo::WireID> > corner::CornerFinderAlg::WireData_IDs
private
std::vector< std::tuple<int,TH2F,int,int> > corner::CornerFinderAlg::WireData_trimmed_histos
private

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