Ellipse2D.ih

/**
 *  This program is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as
 *  published by the Free Software Foundation, either version 3 of the
 *  License, or  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
 *
 **/
 
/**
 * @file Ellipse2D.ih
 * @author David Coeurjolly (\c david.coeurjolly@liris.cnrs.fr )
 * Laboratoire d'InfoRmatique en Image et Systèmes d'information - LIRIS (CNRS, UMR 5205), CNRS, France
 * @author Jacques-Olivier Lachaud (\c jacques-olivier.lachaud@univ-savoie.fr )
 * Laboratory of Mathematics (CNRS, UMR 5807), University of Savoie, France
 *
 * @date 2011/04/12
 *
 * Implementation of inline methods defined in Ellipse2D.h
 *
 * This file is part of the DGtal library.
 */
 
 
//////////////////////////////////////////////////////////////////////////////
#include <cstdlib>
//////////////////////////////////////////////////////////////////////////////
 
///////////////////////////////////////////////////////////////////////////////
// IMPLEMENTATION of inline methods.
///////////////////////////////////////////////////////////////////////////////
 
#define ELLIPSE2D_3_PI_2 (1.5*M_PI)
#define ELLIPSE2D_2_PI   (2.*M_PI)
 
///////////////////////////////////////////////////////////////////////////////
// ----------------------- Standard services ------------------------------
 
template <typename T>
inline
DGtal::Ellipse2D<T>::Ellipse2D(const double x0, const double y0,
             const double a0, const double a1, const double theta)
  : myCenter(x0,y0), myAxis1(a0), myAxis2(a1), myTheta(theta)
{}
 
template <typename T>
inline
DGtal::Ellipse2D<T>::Ellipse2D(const RealPoint &aPoint,
             const double a0, const double a1, const double theta)
  : myCenter(aPoint), myAxis1(a0), myAxis2(a1), myTheta(theta)
{}
 
template <typename T>
inline
DGtal::Ellipse2D<T>::Ellipse2D(const Ellipse2D& other)
  : myCenter(other.myCenter), myAxis1(other.myAxis1),
    myAxis2(other.myAxis2), myTheta(other.myTheta)
{}
 
/////////////////////////////////////////////////////////////////////////////
// ------------- Implementation of 'StarShaped' services ------------------
 
/**
 * @param pp any point in the plane.
 *
 * @return the angle parameter between 0 and 2*Pi corresponding to
 * this point for the shape.
 */
template <typename T>
inline
double
DGtal::Ellipse2D<T>::parameter( const RealPoint & pp ) const
{
  RealPoint v2d( pp - myCenter );
  
  double angle; 
 
  if ( isAlmostEqual(v2d[0],0.) ) 
    angle = v2d[1] > 0. ? M_PI_2 : ELLIPSE2D_3_PI_2;
  else if ( v2d[0] > 0. )
    angle = atan(v2d[1]/v2d[0]) + (v2d[1] < 0. ? ELLIPSE2D_2_PI : 0.);
  else
    angle = atan(v2d[1]/v2d[0]) + M_PI;
 
   return angle;
}
 
/**
 * @param t any angle between 0 and 2*Pi.
 *
 * @return the vector (x(t),y(t)) which is the position on the
 * shape boundary.
 */
template <typename T>
inline
typename DGtal::Ellipse2D<T>::RealPoint 
DGtal::Ellipse2D<T>::x( const double t ) const
{
  const double a2 = myAxis1*myAxis1;
  const double b2 = myAxis2*myAxis2;
  const double costth = cos(t - myTheta);
  const double rho = myAxis2 / sqrt(1. - (a2-b2) / a2 * costth * costth);
 
  return RealPoint(
    rho * cos(t) + myCenter[0],
    rho * sin(t) + myCenter[1] );
}
 
 
/**
 * @param t any angle between 0 and 2*Pi.
 *
 * @return the vector (x'(t),y'(t)) which is the tangent to the
 * shape boundary.
 */
template <typename T>
inline
typename DGtal::Ellipse2D<T>::RealVector 
DGtal::Ellipse2D<T>::xp( const double t ) const
{
  const double a2 = myAxis1 * myAxis1;
  const double b2 = myAxis2 * myAxis2;
  const double costth = cos(t - myTheta);
  const double sintth = sin(t - myTheta);
  const double cost = cos(t);
  const double sint = sin(t);
  const double rho = myAxis2 / sqrt(1. - (a2-b2) / a2 * costth * costth);
  const double rhod = myAxis1 * myAxis2 * (b2-a2) * sintth * costth
    / std::pow(a2 * sintth * sintth + b2 * costth * costth, 1.5);
  
  return RealPoint(
    rhod*cost - rho*sint,
    rhod*sint + rho*cost
  );
}
 
/**
 * @param t any angle between 0 and 2*Pi.
 *
 * @return the vector (x''(t),y''(t)).
 */
template <typename T>
inline
typename DGtal::Ellipse2D<T>::RealVector
DGtal::Ellipse2D<T>::xpp( const double t ) const
{
  const double a2 = myAxis1 * myAxis1;
  const double b2 = myAxis2 * myAxis2;
  const double costth = cos(t - myTheta);
  const double sintth = sin(t - myTheta);
  const double cost = cos(t);
  const double sint = sin(t);
  const double rho = myAxis2 / sqrt( 1. - (a2-b2)/a2*costth*costth);
 
  const double rhod = myAxis1 * myAxis2 * (b2-a2) * sintth * costth
    / std::pow( a2*sintth*sintth + b2*costth*costth, 1.5 );
  const double rhodd = myAxis1 * myAxis2 * (b2-a2)
    / std::pow( a2*sintth*sintth + b2*costth*costth, 2.5 )
    * ( (costth*costth - sintth*sintth) * (a2*sintth*sintth + b2*costth*costth)
        + 3.*(b2-a2)*sintth*sintth*costth*costth );
  
  return RealPoint(
    rhodd*cost - 2.*rhod*sint - rho*cost,
    rhodd*sint + 2.*rhod*cost - rho*sint
  );
}
 
 
///////////////////////////////////////////////////////////////////////////////
// Interface - public :
 
/**
 * Writes/Displays the object on an output stream.
 * @param out the output stream where the object is written.
 */
template <typename T>
inline
void
DGtal::Ellipse2D<T>::selfDisplay ( std::ostream & out ) const
{
  out << "[Ellipse2D] center= " << myCenter
      << " big axis=" << myAxis1
      << " small axis=" << myAxis2
      << " phase=" << myTheta;
}
 
/**
 * Checks the validity/consistency of the object.
 * @return 'true' if the object is valid, 'false' otherwise.
 */
template <typename T>
inline
bool
DGtal::Ellipse2D<T>::isValid() const
{
    return true;
}
 
 
 
///////////////////////////////////////////////////////////////////////////////
// Implementation of inline functions                                        //
 
template <typename T>
inline
std::ostream&
DGtal::operator<< ( std::ostream & out, 
      const Ellipse2D<T> & object )
{
  object.selfDisplay( out );
  return out;
}
 
//                                                                           //
///////////////////////////////////////////////////////////////////////////////