irt/include/IRT.h

0001
0002 #include <vector>
0003
0004 #include <TVector3.h>
0005 #include <TObject.h>
0006 #include <TRef.h>
0007
0008 #ifndef _IRT_FACTORY_
0009 #define _IRT_FACTORY_
0010
0011 // Don't even care about the units here; [mm] in fact; just make it practically small;
0012 #define _IRT_PRECISION_DEFAULT_     (1E-7)
0013 // Well, assume radians in [theta,phi] -> 1urad step should be good enough;
0014 #define _IRT_JACOBIAN_STEP_DEFAULT_ (1E-6)
0015 #define _IRT_ITERATION_LIMIT_        ( 50)
0016 // Here perhaps want something tangible, on a [mm] scale (?), matching the uncertaintly
0017 // in the sensor position measurement, and emission point uncertainty along the trajectory;
0018 #define _IRT_DERIVATIVE_XYZ_STEP_    (1.00)
0019 // Refractive index variation; in units of the (1-n); ~1% makes sense?;
0020 //#define _IRT_DERIVATIVE_NNN_STEP_    (0.01)
0021
0022 #include "CherenkovRadiator.h"
0023 #include "OpticalBoundary.h"
0024 #include "ParametricSurface.h"
0025 #include "IRTSolution.h"
0026
0027 class IRT: public TObject {
0028  public:
0029  IRT(unsigned sector = 0): m_Sector(sector), m_IterationLimit(_IRT_ITERATION_LIMIT_),
0030     m_Precision(_IRT_PRECISION_DEFAULT_), m_JacobianStep(_IRT_JACOBIAN_STEP_DEFAULT_) {};
0031   ~IRT() {};
0032
0033   void AddOpticalBoundary(OpticalBoundary *boundary) {
0034     _m_OpticalBoundaries.push_back(boundary);
0035   };
0036
0037   void SetIterationLimit(unsigned value) { m_IterationLimit = value; };
0038   //void SetPrecision   (double value) { m_Precision = value; };
0039   void SetJacobianStep(double value) { m_JacobianStep = value; };
0040
0041   IRTSolution Solve(const TVector3 &xfrom, const TVector3 &nfrom, const TVector3 &xto, const TVector3 &beam,
0042             //CherenkovRadiator *derivatives = 0, const IRTSolution *seed = 0);
0043             bool derivatives = false, const IRTSolution *seed = 0);
0044   IRTSolution Solve(const TVector3 &xfrom, const TVector3 &nfrom, const double xy[2], const TVector3 &beam,
0045             bool derivatives = false, const IRTSolution *seed = 0);
0046
0047   unsigned GetSector( void ) const { return m_Sector; };
0048
0049   inline OpticalBoundary *tail( void ) const {
0050     return _m_OpticalBoundaries.size() ? GetOpticalBoundary(_m_OpticalBoundaries.size()-1) : 0;
0051   };
0052
0053  private:
0054   bool Transport(const TVector3 &xfrom, const TVector3 &nfrom);
0055   inline OpticalBoundary *GetOpticalBoundary(unsigned id) const {
0056     return (id < _m_OpticalBoundaries.size() ?
0057             dynamic_cast<OpticalBoundary*>(_m_OpticalBoundaries[id].GetObject()) : 0);
0058     //return (id < _m_OpticalBoundaries.size() ? _m_OpticalBoundaries[id] : 0);
0059   };
0060
0061   // FIXME: this is not the right place for this variable;
0062   unsigned m_Sector;
0063
0064   unsigned m_IterationLimit;
0065   double m_Precision, m_JacobianStep;
0066
0067   std::vector<TRef> _m_OpticalBoundaries;
0068   //std::vector<OpticalBoundary *> _m_OpticalBoundaries;
0069
0070 #ifndef DISABLE_ROOT_IO
0071   ClassDef(IRT, 3);
0072 #endif
0073 };
0074
0075 #endif