//-----------------------------------------------------------
// File and Version Information:
// $Id$
//
// Description:
//      Implementation of class PndRiemannHit
//      see PndRiemannHit.hh for details
//
// Environment:
//      Software developed for the PANDA Detector at FAIR.
//
// Author List:
//      Sebastian Neubert    TUM            (original author)
//
//
//-----------------------------------------------------------

// Panda Headers ----------------------

// This Class' Header ------------------
#include "PndRiemannHit.h"
#include "TMath.h"
// C/C++ Headers ----------------------
#include "CbmHit.h"
#include <assert.h>

// Collaborating Class Headers --------
#include "PndRiemannTrack.h"

// Class Member definitions -----------


// scale factor! r=1 => z'=0.5 !!!
// to normalize cm-range 
#define RIEMANNSCALE 1.

ClassImp(PndRiemannHit)


PndRiemannHit::PndRiemannHit(double r, double phi) : fHit(NULL), fCovX(3,3)
{
  double r2=r*r;					//projection to a paraboloid
  fX.SetX(r*TMath::Cos(phi));
  fX.SetY(r*TMath::Sin(phi));
  fX.SetZ(r2);
}

PndRiemannHit::PndRiemannHit(CbmHit* cl) : fHit(cl), fCovX(3,3)
{
  TVector2 a(cl->GetX(),cl->GetY());
  double r=a.Mod()/RIEMANNSCALE; 
  double phi=a.Phi();
  double r2=r*r;
 // fX.SetX(r*TMath::Cos(phi));
 // fX.SetY(r*TMath::Sin(phi));
  fX.SetX(cl->GetX()/RIEMANNSCALE);
  fX.SetY(cl->GetY()/RIEMANNSCALE);
  fX.SetZ(r2);
  
  fCovX[0][0] = TMath::Power(cl->GetDx()/RIEMANNSCALE,2);
  fCovX[1][1] = TMath::Power(cl->GetDy()/RIEMANNSCALE,2);
  fCovX[2][2] = TMath::Power(2*cl->GetX()*cl->GetDx(),2) + TMath::Power(2*cl->GetY()*cl->GetDy(),2)/RIEMANNSCALE;
}


PndRiemannHit::~PndRiemannHit(){}

double
PndRiemannHit::z()const {
  if(fHit==NULL)return 0;
  else return fHit->GetZ();
}

void 
PndRiemannHit::calcPosOnTrk(PndRiemannTrack* trk){
  assert(trk!=NULL);
  TVectorD o=trk->orig(); 
  o*=RIEMANNSCALE;// convert back to cm;
  double r=trk->r();
  r*=RIEMANNSCALE;
  const CbmHit* firstHit=trk->getHit(0)->hit();
  assert(firstHit!=NULL);
  TVector2 k(firstHit->GetX()-o[0],firstHit->GetY()-o[1]);
  TVector2 l(fHit->GetX()-o[0],fHit->GetY()-o[1]);
  fAlpha=l.DeltaPhi(k);
  fS=fAlpha*r;
}

double PndRiemannHit::sigmaXY() const
{
	double result = TMath::Sqrt(fCovX[0][0] + fCovX[1][1]);
	if (result == 0.){
		return 1.;
	}
	else return result;
}