#ifndef PAMVMCPRIMARYGF
#define PAMVMCPRIMARYGF

#include <vector>
//#include <TObjectTable.h>
#include "TObject.h"
#include "TrkParams.h"
#include "TMath.h"
#define ZOFFSET 49.229


using TMath::Min;

struct PamVMCPlCross
{
  Double_t fX, fY;
  Int_t fplID;
  PamVMCPlCross(Double_t x,Double_t y, Int_t id):fX(x), fY(y), fplID(id){ };
  virtual ~PamVMCPlCross(){; };
};

typedef vector<PamVMCPlCross*> Plcrvec;


//variable NGF is defined in TrkParams.h
class  PamVMCPrimaryGF : public TObject
{
  Plcrvec vec;
  Double_t fXtol, fYtol;
  Bool_t fisinside;
  
  Int_t PlCross(Double_t z_prev, Double_t z_curr){
    for(Int_t i = 0; i<NGF; i++){
      if( ( (z_prev>=TrkParams::zGF[i] + ZOFFSET ) && (z_curr<=TrkParams::zGF[i] +  ZOFFSET ) ) || ( (z_prev<=TrkParams::zGF[i] + ZOFFSET ) && (z_curr>=TrkParams::zGF[i] +  ZOFFSET ) ) ) {
      return i;

    }
  }
  return -1;
}
  
 Bool_t IsInside(Int_t id, Double_t z_prev, Double_t x_prev, Double_t y_prev, Double_t z_curr, Double_t x_curr, Double_t y_curr, Double_t* xc, Double_t* yc){
   if (z_prev==z_curr){
     *xc = x_curr;
     *yc = y_curr;
     if( (TrkParams::xGF_min[id]<=x_curr) && (TrkParams::xGF_max[id]>=x_curr) && 
	 (TrkParams::yGF_min[id]<=y_curr) && (TrkParams::yGF_max[id]>=y_curr) ) return kTRUE; 
   } else {
     Double_t kx = (x_prev - x_curr)/(z_prev-z_curr);
     Double_t bx = x_prev - kx*z_prev;
     Double_t ky = (y_prev - y_curr)/(z_prev-z_curr);
     Double_t by = y_prev - ky*z_prev;
     *xc = kx*(TrkParams::zGF[id]+ZOFFSET)+bx;
     *yc = ky*(TrkParams::zGF[id]+ZOFFSET)+by;
     //cout<<"Trk_plane:"<<id<<"  x_prev,y_prev,zprev:"<<x_prev<<","<<y_prev<<","<<z_prev<<endl;
     //cout<<"Trk_plane:"<<id<<"  x_curr,y_curr,zcurr:"<<x_curr<<","<<y_curr<<","<<z_curr<<endl;
     //cout<<"Trk_plane:"<<id<<"  xc,yc,zc:"<<*xc<<","<<*yc<<","<<TrkParams::zGF[id]+ZOFFSET<<endl;
     if( (TrkParams::xGF_min[id]<=*xc) && (TrkParams::xGF_max[id]>=*xc) && 
	 (TrkParams::yGF_min[id]<=*yc) && (TrkParams::yGF_max[id]>=*yc) ) return kTRUE; 
   }
   return kFALSE;
 }

 void AddPlane(Double_t x_cr, Double_t y_cr, Int_t pl_id){
   vec.push_back( new PamVMCPlCross( x_cr, y_cr, pl_id) );
 }
 


 public:
 PamVMCPrimaryGF():fXtol(-100.),fYtol(-100.), fisinside(kFALSE){  };
 virtual ~PamVMCPrimaryGF() {  };


 //doing after each step of primary
 void CheckCrossing(Double_t z_prev, Double_t x_prev, Double_t y_prev, Double_t z_curr, Double_t x_curr, Double_t y_curr){
   Int_t pl_id = PlCross(z_prev, z_curr);
   Double_t xc,yc;
   xc = yc = -100.;
   if(pl_id >= 0){
     fisinside = IsInside(pl_id, z_prev, x_prev, y_prev, z_curr,  x_curr, y_curr, &xc, &yc);
     if (fisinside) AddPlane(xc, yc, pl_id);
   }
 }
 //doing after finishing primary track
 Bool_t IsInsideAcceptance(Double_t* xc, Double_t* yc, Double_t &xtol, Double_t &ytol){
   Plcrvec::const_iterator p = vec.begin();
   *xc = *yc = -100;
   xtol = ytol = 100.;
   Bool_t is_tol_outside = kFALSE;
    while( p!=vec.end() )
      {
        xc[(*p)->fplID] = (*p)->fX;
	yc[(*p)->fplID] = (*p)->fY;
	Double_t xtolpl = Min((*p)->fX-TrkParams::xGF_min[(*p)->fplID], TrkParams::xGF_max[(*p)->fplID]-(*p)->fX);
	Double_t ytolpl = Min((*p)->fY-TrkParams::yGF_min[(*p)->fplID], TrkParams::yGF_max[(*p)->fplID]-(*p)->fY);
	if( (xtolpl < 0.) || (ytolpl<0.)) is_tol_outside = kTRUE;
	//cout<<"xc="<<(*p)->fX<<" yc="<<(*p)->fY<<" xtolpl="<<xtolpl<<" ytolpl="<<ytolpl<<endl;
	xtol=(xtolpl<xtol)?xtolpl:xtol;
	ytol=(ytolpl<ytol)?ytolpl:ytol;
        p++;
      }
    if ( is_tol_outside ) xtol = ytol = -100.;
    //cout<<"xtolpl_FIN="<<xtol<<" ytolpl_FIN="<<ytol<<" PL:"<<vec.size()<<endl;
    if (vec.size() == NGF) return kTRUE;
    xtol = ytol = -100.;
    return kFALSE;
 }

 // calls after primary propagation to clean everything
 virtual void Clear(Option_t * = ""){
   fXtol=fYtol=-100.;
   vec.clear();
 }
};
#endif
