#include "PhysicsVectors/HepThreeVector.h"
#include "PhysicsVectors/HepRotation.h"
#include "PhysicsVectors/HepLorentzVector.h"
#include "PhysicsVectors/HepLorentzRotation.h"

#include "ZMutility/iostream"
using std::cout;
using std::cerr;
using std::endl;

ZM_USING_NAMESPACE( zmpv )


void checkCLHEProtate (double psi, UnitVector axis)  {
  cout << "axis = " << axis << " psi = " << psi << endl;

  double theta = axis.theta();
  double phi = axis.phi();
  Rotation R ( axis, psi );
  cout << "R = " << R;

  Rotation Rclhep = RotationZ(phi)*RotationY(theta)
                *RotationZ(psi)*RotationY(-theta)*RotationZ(-phi);
  cout << "Rclhep = " << Rclhep;

  if ( R.isNear(Rclhep) )
    cout << "OK\n";
  else
    cout << "Discrepancy:\nRclhep is really " << Rclhep.axisAngle() << endl;
}


int  main()  {

  UnitVector axis;
  axis.set (2, 4, 6);
  checkCLHEProtate ( .5, axis);
  checkCLHEProtate ( -.23, axis);
  axis.set (-2, 6, 5);
  checkCLHEProtate ( 2.5, axis);
  checkCLHEProtate ( -2.8, axis);

  return  0;

}  // main()