JRECONSTRUCTION::JShowerDirectionPrefit Class Reference

class to handle the direction fit of the shower reconstruction, mainly dedicated for ORCA More...

#include <JShowerDirectionPrefit.hh>

Inheritance diagram for JRECONSTRUCTION::JShowerDirectionPrefit:

Classes
struct	input_type
	Input data type. More...

Public Types
typedef JRegressor< JShower3EZ, JAbstractMinimiser >	JRegressor_t
typedef JModuleL0	module_type
typedef std::vector< module_type >	detector_type

Public Member Functions
std::vector< double >	computeChi2 (const std::vector< double > &Ev, const std::vector< JFIT::JPMTW0 > &data, const JGEOMETRY3D::JRotation3D &R, const JShower3EZ &sh)
	Compute chi2 for a set of energies.
std::pair< double, double >	findBestEnergy (std::vector< double > Ev, const std::vector< JFIT::JPMTW0 > &data, const JGEOMETRY3D::JRotation3D &R, const JShower3EZ &sh)
	Find best energy iteratively increasing the upper threshold if there is suspicion the best-fit is outside.
	JShowerDirectionPrefit (const JShowerDirectionPrefitParameters_t &parameters, const storage_type &storage, const JPMTParametersMap &pmtParameters, const int debug=0)
	Parameterized constructor.
input_type	getInput (const JModuleRouter &router, const JSummaryRouter &summary, const JDAQEvent &event, const JEvt &in, const coverage_type &coverage) const
	Get input data.
JEvt	operator() (const input_type &input)
	Fit function.
void	reset ()
	Reset fit parameters.
bool	equals (const JShowerDirectionPrefitParameters_t &parameters) const
	Equality.
	ClassDef (JShowerDirectionPrefitParameters_t, 2)

Public Attributes
std::vector< double >	Ev
const JPMTParametersMap &	pmtParameters
size_t	numberOfPrefits
	number of prefits
double	TMax_ns
	maximum time for local coincidences [ns]
double	TMin_ns
	minimum time for local coincidences [ns]
double	DMax_m
	maximal distance to optical module [m]
double	R_Hz
	default rate [Hz]
double	VMax_npe
	maximum number of of photo-electrons
double	scanAngle_deg
	scanning angle step in [deg]
double	Emin_GeV
	minimum energy to scan
double	Emax_GeV
	maximum energy to scan
int	En
	number of points to scan in energy range

Detailed Description

class to handle the direction fit of the shower reconstruction, mainly dedicated for ORCA

Definition at line 70 of file JShowerDirectionPrefit.hh.

Member Typedef Documentation

JRegressor_t

JRegressor<JShower3EZ, JAbstractMinimiser> JRECONSTRUCTION::JShowerDirectionPrefit::JRegressor_t

Definition at line 75 of file JShowerDirectionPrefit.hh.

module_type

JModuleL0 JRECONSTRUCTION::JShowerDirectionPrefit::module_type

Definition at line 76 of file JShowerDirectionPrefit.hh.

detector_type

std::vector<module_type> JRECONSTRUCTION::JShowerDirectionPrefit::detector_type

Definition at line 77 of file JShowerDirectionPrefit.hh.

Constructor & Destructor Documentation

JShowerDirectionPrefit()

JRECONSTRUCTION::JShowerDirectionPrefit::JShowerDirectionPrefit ( const JShowerDirectionPrefitParameters_t & parameters, const storage_type & storage, const JPMTParametersMap & pmtParameters, const int debug = 0 )

inline

Parameterized constructor.

Parameters

parameters	parameters
storage	storage
pmtParameters	PMT parameters
debug	debug

Definition at line 222 of file JShowerDirectionPrefit.hh.

                                                                                 :
      JShowerDirectionPrefitParameters_t(parameters),
      JRegressor_t(storage),
      pmtParameters(pmtParameters)
    {
      using namespace JPP;
      
      JRegressor_t::debug  = debug;
      JRegressor_t::T_ns.setRange(parameters.TMin_ns, parameters.TMax_ns);
      JRegressor_t::Vmax_npe = VMax_npe;
          
      if (Emin_GeV > Emax_GeV || En <= 1) {
        THROW(JException, "Invalid energy input " << Emin_GeV << ' ' << Emax_GeV << ' ' << En);
      }
 
      const double base = std::pow((Emax_GeV / Emin_GeV), 1.0 / (En - 1));
      
      for (int i = 0; i != En; ++i) {
        Ev.push_back(Emin_GeV * std::pow(base, i));
      }
    }

Member Function Documentation

computeChi2()

std::vector< double > JRECONSTRUCTION::JShowerDirectionPrefit::computeChi2 ( const std::vector< double > & Ev, const std::vector< JFIT::JPMTW0 > & data, const JGEOMETRY3D::JRotation3D & R, const JShower3EZ & sh )

inline

Compute chi2 for a set of energies.

Parameters

Ev	energy vector
data	PMTs
R	Rotation
sh	shower

Definition at line 126 of file JShowerDirectionPrefit.hh.

                                      {
      using namespace JPP;
 
      std::vector<double> chi2(Ev.size(), 0.0);
 
      for (auto i = data.begin(); i != data.end(); ++i) {
 
        JPMTW0 pmt = *i;
        pmt.rotate(R);
 
        auto H1 = getH1(sh, pmt);
        auto H0 = getH0(pmt.getR());
 
        const bool hit = pmt.getN() != 0;
        const double chi2_bg = getChi2(get_value(H0), hit);
 
        for (size_t j = 0; j < Ev.size(); ++j) {
          chi2[j] += getChi2(get_value(Ev[j] * H1  +  H0), hit) - chi2_bg; ; // H1 is linear with E, -log-lik ratio
        }
      }
 
      return chi2;
    }

findBestEnergy()

std::pair< double, double > JRECONSTRUCTION::JShowerDirectionPrefit::findBestEnergy ( std::vector< double > Ev, const std::vector< JFIT::JPMTW0 > & data, const JGEOMETRY3D::JRotation3D & R, const JShower3EZ & sh )

inline

Find best energy iteratively increasing the upper threshold if there is suspicion the best-fit is outside.

Parameters

Ev	energy vector
data	PMTs
R	Rotation
sh	shower

Definition at line 163 of file JShowerDirectionPrefit.hh.

                                               {
      using namespace JPP;
 
      const int max_extensions = 5; 
      int extension_count = 0;
 
      while (true) {
 
        std::vector<double> chi2 = computeChi2(Ev, data, R, sh);
 
        auto p = std::min_element(chi2.begin(), chi2.end());
        size_t idx = p - chi2.begin();
 
        bool at_upper_edge = (idx == Ev.size() - 1);
 
        bool monotonic_tail = false;
        if (Ev.size() >= 3) {
          size_t n = Ev.size();
          monotonic_tail = (chi2[n-3] > chi2[n-2]) &&
            (chi2[n-2] > chi2[n-1]);
        }
 
        if (!(at_upper_edge && monotonic_tail)) {
          return {Ev[idx], *p};
        }
 
        if (extension_count >= max_extensions) {
          return {Ev[idx], *p};
        }
 
        extension_count++;
 
        std::vector<double> Ev_ext;
        double E_start = Ev.back();
        double E_end   = 100.0 * E_start;
 
        size_t N = Ev.size();
 
        for (size_t k = 0; k < N; ++k) {
          double t = double (k)  / (N - 1);
          Ev_ext.push_back(E_start * std::pow(E_end / E_start, t));
        }
 
        Ev = std::move(Ev_ext);  
      }
    }

getInput()

input_type JRECONSTRUCTION::JShowerDirectionPrefit::getInput ( const JModuleRouter & router, const JSummaryRouter & summary, const JDAQEvent & event, const JEvt & in, const coverage_type & coverage ) const

inline

Get input data.

Parameters

router	module router
summary	summary data
event	event
in	start values
coverage	coverage

Returns

input data

Definition at line 257 of file JShowerDirectionPrefit.hh.

    {
      using namespace std;
      using namespace JTRIGGER;
 
      input_type   input(event.getDAQEventHeader(), in, coverage);
 
      const JBuildL0 <JHitR0>   buildL0;
      map<int, vector<JHitR0> > data;
 
      const JDAQTimeslice timeslice(event, true);
 
      JSuperFrame2D<JHit> buffer;
 
      for (JDAQTimeslice::const_iterator i = timeslice.begin(); i != timeslice.end(); ++i) {
 
        if (router.hasModule(i->getModuleID())) {
 
          buffer(*i, router.getModule(i->getModuleID()));
 
          buildL0(buffer, back_inserter(data[i->getModuleID()]));
        }
      }
 
      for (const auto& module : router.getReference()) {
        if (!module.empty()) {
          input.data.push_back(module_type(module, summary.getSummaryFrame(module.getID(), R_Hz), data[module.getID()]));
        }
      }
 
      return input;
    }

operator()()

JEvt JRECONSTRUCTION::JShowerDirectionPrefit::operator() ( const input_type & input )

inline

Fit function.

Parameters

input

input data

Returns

fit results

Definition at line 300 of file JShowerDirectionPrefit.hh.

    {
      using namespace std;
      using namespace JPP;
 
      JEvent event(JSHOWERDIRECTIONPREFIT);
 
      JEvt out;
 
      JEvt in = input.in;
      
      in.select(numberOfPrefits, qualitySorter);
 
      if (!in.empty()) {
        in.select(JHistory::is_event(in.begin()->getHistory()));
      }
  
      for (JEvt::const_iterator shower = in.begin(); shower != in.end(); ++shower) {
 
        vector<JPMTW0> data;
 
        const JPosition3D vertex(getPosition(*shower)); 
        const double      time = shower->getT(); 
        
        for (const auto& module : input.data) {
          
          JPosition3D pos(module->getPosition());
 
          pos.sub(vertex);
          
          if (pos.getLength() <= DMax_m) {
            
            const double t1 = time + pos.getLength() * getInverseSpeedOfLight() * getIndexOfRefraction();
 
            for (size_t i = 0; i != module->size(); ++i) {
 
              if (module.getStatus(i)) {
 
                struct {
 
                  bool operator()(const JHitR0& hit) const
                  {
                    return (hit.getPMT() == pmt && T_ns(hit.getT()));
                  }
                    
                  const JTimeRange T_ns;
                  const size_t     pmt;
 
                } match = { JRegressor_t::T_ns + t1, i };
 
                const JPMTIdentifier id(module->getID(), i);
 
                const JPMTParameters& wip = pmtParameters.getPMTParameters(id);
 
                const size_t ns = count_if(module.begin(), module.end(), match);
                const double QE = wip.QE;
 
                JPMT pmt = module->getPMT(i);
 
                pmt.sub(vertex);
 
                data.push_back(JPMTW0(pmt, QE, module.frame.getRate(i), ns));
                
              }
            }
          }
        }
 
        JVector3D start_dir(0,0,0);
 
        for (vector<JPMTW0>::const_iterator i = data.begin(); i != data.end(); ++i) {
          if (i->getN() != 0) start_dir.add(i->getPosition());
        }
        
        JOmega3D scan_directions(JDirection3D(start_dir),
                                 JOmega3D::range_type(0.0,PI),
                                 scanAngle_deg * PI / 180.0);
        
        const JShower3EZ sh(JVertex3D(JVector3D(0,0,0), shower->getT()), JVersor3Z(), 1);
 
        for (JOmega3D_t::const_iterator dir = scan_directions.begin(); dir != scan_directions.end(); ++dir) {
 
          const JGEOMETRY3D::JRotation3D R(*dir);
          
          std::pair<double,double> result = findBestEnergy(Ev, data, R, sh);
          
          out.push_back(getFit(JHistory(shower->getHistory(), event()),
                               JShower3D(JVertex3D(vertex,time), JDirection3D(*dir)),
                               getQuality(result.second),
                               data.size(),
                               result.first));
 
          // set additional values
          out.rbegin()->setW(JPP_COVERAGE_ORIENTATION, input.coverage.orientation);
          out.rbegin()->setW(JPP_COVERAGE_POSITION,    input.coverage.position);
 
        }
      }
      
      // apply default sorter
 
      sort(out.begin(), out.end(), qualitySorter);
 
      copy(input.in.begin(), input.in.end(), back_inserter(out));
 
      return out;
    
    }

reset()

void JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::reset ( )

inlineinherited

Reset fit parameters.

Definition at line 35 of file JShowerDirectionPrefitParameters_t.hh.

    {
      numberOfPrefits             =    1;
      TMax_ns                     =   30;
      TMin_ns                     =  -30;
      DMax_m                      =   80;
      scanAngle_deg               =   10;
      R_Hz                        = 10.0e3;
      VMax_npe                    =   20.0;
      Emin_GeV                    = 1;
      Emax_GeV                    = 1000;
      En                          = 50;
    }

equals()

bool JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::equals ( const JShowerDirectionPrefitParameters_t & parameters ) const

inlineinherited

Equality.

Parameters

parameters

fit parameters

Returns

true if equals; else false

Definition at line 55 of file JShowerDirectionPrefitParameters_t.hh.

    {
      return (this->TMax_ns            == parameters.TMax_ns                  &&         
              this->TMin_ns                  == parameters.TMin_ns                  &&         
              this->numberOfPrefits          == parameters.numberOfPrefits          &&                 
              this->DMax_m                   == parameters.DMax_m                   &&
              this->R_Hz                     == parameters.R_Hz                     && 
              this->scanAngle_deg            == parameters.scanAngle_deg            &&
              this->VMax_npe                 == parameters.VMax_npe                 &&
              this->Emin_GeV                 == parameters.Emin_GeV                 &&
              this->Emax_GeV                 == parameters.Emax_GeV                 &&
              this->En                       == parameters.En                         ); 
    }

ClassDef()

JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::ClassDef ( JShowerDirectionPrefitParameters_t , 2  )

inherited

Member Data Documentation

Ev

std::vector<double> JRECONSTRUCTION::JShowerDirectionPrefit::Ev

Definition at line 83 of file JShowerDirectionPrefit.hh.

pmtParameters

const JPMTParametersMap& JRECONSTRUCTION::JShowerDirectionPrefit::pmtParameters

Definition at line 409 of file JShowerDirectionPrefit.hh.

numberOfPrefits

size_t JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::numberOfPrefits

inherited

number of prefits

Definition at line 71 of file JShowerDirectionPrefitParameters_t.hh.

TMax_ns

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::TMax_ns

inherited

maximum time for local coincidences [ns]

Definition at line 72 of file JShowerDirectionPrefitParameters_t.hh.

TMin_ns

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::TMin_ns

inherited

minimum time for local coincidences [ns]

Definition at line 73 of file JShowerDirectionPrefitParameters_t.hh.

DMax_m

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::DMax_m

inherited

maximal distance to optical module [m]

Definition at line 74 of file JShowerDirectionPrefitParameters_t.hh.

R_Hz

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::R_Hz

inherited

default rate [Hz]

Definition at line 75 of file JShowerDirectionPrefitParameters_t.hh.

VMax_npe

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::VMax_npe

inherited

maximum number of of photo-electrons

Definition at line 76 of file JShowerDirectionPrefitParameters_t.hh.

scanAngle_deg

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::scanAngle_deg

inherited

scanning angle step in [deg]

Definition at line 77 of file JShowerDirectionPrefitParameters_t.hh.

Emin_GeV

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::Emin_GeV

inherited

minimum energy to scan

Definition at line 78 of file JShowerDirectionPrefitParameters_t.hh.

Emax_GeV

double JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::Emax_GeV

inherited

maximum energy to scan

Definition at line 79 of file JShowerDirectionPrefitParameters_t.hh.

En

int JRECONSTRUCTION::JShowerDirectionPrefitParameters_t::En

inherited

number of points to scan in energy range

Definition at line 80 of file JShowerDirectionPrefitParameters_t.hh.

---

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

• JShowerDirectionPrefit.hh