#ifndef __JGEOMETRY3DTESTKIT__
#define __JGEOMETRY3DTESTKIT__

#include "TRandom3.h"

#include "JGeometry3D/JPosition3D.hh"
#include "JGeometry3D/JDirection3D.hh"
#include "JGeometry3D/JRotation3D.hh"
#include "JGeometry3D/JQuaternion3D.hh"
#include "JGeometry3D/JEulerMatrix3D.hh"
#include "JGeometry3D/JSphere3D.hh"
#include "JGeometry3D/JCylinder3D.hh"
#include "JMath/JRandom.hh"
#include "JMath/JConstants.hh"

/**
 * \author mdejong
 */

namespace JGEOMETRY3D {}
namespace JPP { using namespace JGEOMETRY3D; }

namespace JGEOMETRY3D {

  using JMATH::PI;
  using JMATH::getRandom;


  /**
   * Randomize position.
   *
   * \param  p             pointer to valid object
   */
  inline void randomize(JPosition3D* p)
  {
    new (p) JPosition3D(getRandom<double>(-1.0, +1.0),
			getRandom<double>(-1.0, +1.0),
			getRandom<double>(-1.0, +1.0));
  };


  /**
   * Randomize direction.
   *
   * \param  p             pointer to valid object
   */
  inline void randomize(JDirection3D* p)
  {
    double x, y, z;

    gRandom->Sphere(x, y, z, 1.0);

    new (p) JDirection3D(x, y, z);
  };


  /**
   * Randomize quaternion.
   *
   * \param  p             pointer to valid object
   */
  inline void randomize(JQuaternion3D* p)
  {
    double x, y, z;

    gRandom->Sphere(x, y, z, 1.0);

    const double alpha = getRandom<double>(-2.0*PI, +2.0*PI);

    new (p) JQuaternion3D(alpha, JVersor3D(x, y, z));
  }


  /**
   * Randomize Euler angle.
   *
   * \param  p             pointer to valid object
   */
  inline void randomize(JEulerAngle3D* p)
  {
    new (p) JEulerAngle3D(getRandom<double>(0.0, PI),
			  getRandom<double>(0.0, PI),
			  getRandom<double>(0.0, PI));
  }


  /**
   * Randomize Euler matrix.
   *
   * \param  p             pointer to valid object
   */
  inline void randomize(JEulerMatrix3D* p)
  {
    new (p) JEulerMatrix3D(getRandom<JEulerAngle3D>());
  }


  /**
   * Get random position.\n
   * The positions are uniformly distributed in given sphere.
   *
   * \param  sphere       sphere
   * \return              position
   */ 
  inline JVector3D getRandomPosition(const JSphere3D& sphere)
  {
    const double R3max = sphere.getRadius() * sphere.getRadius() * sphere.getRadius();

    const double R3  = gRandom->Uniform(0.0, R3max);
    const double ct  = gRandom->Uniform(-1.0, +1.0);
    const double phi = gRandom->Uniform(-PI, +PI);
    
    const double R   = cbrt(R3);
    const double st  = sqrt((1.0 + ct) * (1.0 - ct));

    return sphere.getPosition() + JVector3D(R * st * cos(phi),
					    R * st * sin(phi),
					    R * ct);
  }


  /**
   * Get random position.\n
   * The positions are uniformly distributed in given cylinder.
   *
   * \param  cylinder     cylinder
   * \return              position
   */ 
  inline JVector3D getRandomPosition(const JCylinder3D& cylinder)
  {
    const double R2max = cylinder.getRadius() * cylinder.getRadius();

    const double R2  = gRandom->Uniform(0.0, R2max);
    const double phi = gRandom->Uniform(-PI, +PI);
    const double z   = gRandom->Uniform(cylinder.getZmin(), cylinder.getZmax());

    const double R   = sqrt(R2);

    return JVector3D(cylinder.getX()  +  R * cos(phi),
		     cylinder.getY()  +  R * sin(phi),
		     z);
  }
}

#endif