//*-- Authors: A.Ivashkin && A.Sadovsky (04.11.2004)
#ifndef HPARTICLERK_H
#define HPARTICLERK_H

#include "hmdctrackgfield.h"

#include "TCollection.h"     // kIterForward + kIterBackward
#include "TObjArray.h"
#include "TString.h"
#include "TObject.h"
#include "TVector3.h"
#include "TVectorD.h"




class HParticleRKPlane : public TObject {

 private:
    TVector3 vCenter; //! A point on the plane.
    TVector3 vNormal; //! Normal vector of the plane.

 public:

    HParticleRKPlane() {  }
    HParticleRKPlane(const TVector3 &center, const TVector3 &normal);

    virtual ~HParticleRKPlane();

    virtual Double_t distanceToPlane       (const TVector3 &point) const;
    virtual Bool_t   findIntersection      (TVector3 &pointIntersect, const TVector3 &pos, const TVector3 &dir) const;
    virtual Bool_t   isOnSurface           (const TVector3 &point) const;
    virtual void     print                 (Option_t* opt="") const;
    virtual Double_t signedDistanceToPlane (const TVector3 &point) const;
    virtual void     transform             (const TRotation &transMat);

    virtual TVector3 const& getCenter() const { return vCenter; }
    virtual TVector3 const& getNormal() const { return vNormal; }
    virtual Bool_t          setPlane (const TVector3 &origin, const TVector3 &normal);

    ClassDef(HParticleRKPlane, 0)
};




//----------------------------------------------------------------------------








class HParticleRK : public TObject {
 private:
   HMdcTrackGField* fieldMap;              // pointer to field map
   Double_t  fieldScaleFact;               //! Factor of the magnetic field strength.

   //--------------------------------- RK-tracking ----------------------------
   Double_t   initialStepSize;               // initial RK step size
   Double_t   stepSizeInc;                   // factor to increase stepsize
   Double_t   stepSizeDec;                   // factor to decrease stepsize
   Double_t   maxStepSize;                   // maximum stepsize
   Double_t   minStepSize;                   // minimum stepsize
   Double_t   minPrecision;                  // minimum required precision for a single step
   Double_t   maxPrecision;                  // maximum required precision for a single step
   Int_t      maxNumSteps;                   // maximum number of steps
   Double_t   maxDist;                       // maximum distance for straight line approximation

   Double_t   trackLength;                   // track length

   Int_t      jstep;                         // step number

   Bool_t     bFillPointArrays; // == kTRUE if point arrays should be filled
   Bool_t     bOwnArrays;       // == kTRUE delete Arrays each progagate
   TObjArray  pointsTrack;      // Points along trajectory stored
   TObjArray  fieldTrack;       // Field at Points along trajectory stored

   void     calcField_mm(const TVector3& xv, TVector3& btos, Double_t fpol) const;
   Double_t rkgtrk(Double_t,Double_t,TVector3&,TVector3&);
   void     rkeqfw(TVector3&,Double_t,TVector3&,TVector3&);

public:
   HParticleRK();
   virtual ~HParticleRK() {}
   void             setFieldFact(Double_t fpol)                           { fieldScaleFact = fpol; }
   void             setFieldMap (HMdcTrackGField *fpField, Double_t fpol) { fieldMap = fpField; fieldScaleFact = fpol; }
   void             setFillPointsArrays(Bool_t fill,Bool_t owner)         { bFillPointArrays = fill; bOwnArrays = owner;}
   void             clearPointsArrays()                                   { pointsTrack.Delete(); fieldTrack.Delete();}
   TObjArray const& getPointsTrack()    const                             { return pointsTrack; }
   TVector3         calcFieldIntegral() const;
   Double_t         getTrackLength()                                      { return trackLength; } // RK-track length
   Bool_t           propagateToPlane(const Int_t chrg,TVectorD &stateTo, const TVectorD &stateFrom,const HParticleRKPlane &planeTo,Double_t& length, Bool_t propDir = kIterForward);
   Bool_t           propagateToPCA  (const Int_t chrg,TVectorD &stateTo, const TVectorD &stateFrom,const TVector3& pointTo,Double_t& length, Bool_t propDir = kIterBackward);
   ClassDef(HParticleRK,0)
};


#endif