---

guge_linear.cpp

Go to the documentation of this file.

/* LIBGUL - Geometry Utility Library
 * Copyright (C) 1998-1999 Norbert Irmer
 *
 * This library is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Library General Public
 * License as published by the Free Software Foundation; either
 * version 2 of the License, or (at your option) any later version.
 *
 * This library 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
 * Library General Public License for more details.
 *
 * You should have received a copy of the GNU Library General Public
 * License along with this library; if not, write to the
 * Free Software Foundation, Inc., 59 Temple Place - Suite 330,
 * Boston, MA 02111-1307, USA.
 */

#include "stdafx.h"

#include "gul_types.h"
#include "gul_float.h"
#include "gul_vector.h"
#include "guge_normalize.h"


namespace guge {

using gul::Ptr;
using gul::rtr;
using gul::hpoint;
using gul::point;
using gul::point2;
using gul::normalize;
using gul::cross_product;
using gul::abs_equal;
using gul::rel_equal;
using gul::zero_sine;
using gul::euclid;

/*--------------------------------------------------------------------------
  check, if a curve (given by its homogeneous control points), can be 
  approximated with a line through its endpoints
---------------------------------------------------------------------------*/
template< class T, class HP, class EP >
bool IsLinear( const int n, const Ptr< HP >& Pw, 
               const T tol )
{
  EP A,B,b,u0,X;
  point<T> L,H,R;
  int i;
  bool bzero;

  A = euclid( Pw[0] );
  B = euclid( Pw[n] );

  b = B - A;
  bzero = rel_equal(B,A,rtr<T>::zero_tol());
  
  if( !bzero )
  {
    normalize( &u0, b );
    L = cross_product( u0, A );

    for( i = 1; i < n; i++ )
    {
      X = euclid( Pw[i] );

      H = cross_product( u0, X );
      R = H - L;

      if( !abs_equal( H, L, tol ) )
        return false;
    }
  }
  else
  {
    for( i = 0; i <= n; i++ )
    {
      X = euclid( Pw[i] );
      B = X - A;

      if( !abs_equal(X,A,tol) )
        return false;
    }
  }
  return true;
}
// template instantiation
template bool IsLinear<float,hpoint<float>,point<float> >(
                        const int n, const Ptr< hpoint<float> >& Pw, 
                        const float tol );
template bool IsLinear<double,hpoint<double>,point<double> >(
                        const int n, const Ptr< hpoint<double> >& Pw, 
                        const double tol );
template bool IsLinear<float,point<float>,point<float> >(
                        const int n, const Ptr< point<float> >& Pw, 
                        const float tol );
template bool IsLinear<double,point<double>,point<double> >(
                        const int n, const Ptr< point<double> >& Pw, 
                        const double tol );

template bool IsLinear<float,hpoint2<float>,point2<float> >(
                        const int n, const Ptr< hpoint2<float> >& Pw, 
                        const float tol );
template bool IsLinear<double,hpoint2<double>,point2<double> >(
                        const int n, const Ptr< hpoint2<double> >& Pw, 
                        const double tol );
template bool IsLinear<float,point2<float>,point2<float> >(
                        const int n, const Ptr< point2<float> >& Pw, 
                        const float tol );
template bool IsLinear<double,point2<double>,point2<double> >(
                        const int n, const Ptr< point2<double> >& Pw, 
                        const double tol );


/*-------------------------------------------------------------------------
  check, if a surface can be approximated with two triangles through its
  four corners
---------------------------------------------------------------------------*/
template< class T, class HP >
bool IsPlanar( const int nu, const int nv, const Ptr< Ptr< HP > >& Pw,
               const T tol )
{
  point<T> A,B,C,D,b,c,d,h,n0,u0,X,H,L,R;
  T r,l;
  bool is_point,is_plane, bzero,czero,dzero;
  int i,j;
  
  A = euclid( Pw[0][0] );
  B = euclid( Pw[nv][0] );
  C = euclid( Pw[nv][nu] );
  D = euclid( Pw[0][nu] );

  b = B - A;
  bzero = rel_equal( B, A, rtr<T>::zero_tol() );
  c = C - A;
  czero = rel_equal( C, A, rtr<T>::zero_tol() );
  d = D - A;
  dzero = rel_equal( D, A, rtr<T>::zero_tol() );      

  is_point = false;
  is_plane = false;
  
  if( bzero && czero && dzero )               /* Punkt */
    is_point = true;
  else                                                             /* Ebene */
  {
    if( (!bzero) && (!czero) )
    {
      h = cross_product( b, c );
      is_plane = !zero_sine( b, c, rtr<T>::zero_tol() );  /* b || c ??? */      
    }
       
    if( is_plane )
    {
      normalize( &n0, h );  
    }  
    else
    {  
      if( (!bzero) && (!dzero) )
      {
        h = cross_product( b, d );
        is_plane = !zero_sine( b, d, rtr<T>::zero_tol() );  /* b || d ??? */ 
      }               

      if( is_plane )
      {
        normalize( &n0, h );
      }
      else
      {
        if( (!czero) && (!dzero) )
        {
          h = cross_product( c, d );
          is_plane = !zero_sine( c, d, rtr<T>::zero_tol() );  /* c || d ??? */ 
        }               

        if( is_plane )
        {
          normalize( &n0, h );
        }
        else                                                     /* Gerade */ 
          if( !bzero )
            normalize( &u0, b );
          else
            if( !czero )
              normalize( &u0, c );
            else
              normalize( &u0, d );              
      }
    }            
  } 

  if( is_plane )
  {
    r = n0 * A;
    if( r < 0 )
    {
      r = -r;
      n0 = -n0;
    }  
    
    for( j = 0; j <= nv; j++ )
      for( i = 0; i <= nu; i++ )
      {
        X = euclid( Pw[j][i] );

        l = n0 * X;
        if( gul::rtr<T>::fabs( l - r ) > tol ) return false;
      }
  }                 
  else if( !is_point )
  {
    L = cross_product( u0, A );

    for( j = 0; j <= nv; j++ )
      for( i = 0; i <= nu; i++ )
      {
        X = euclid( Pw[j][i] );

        H = cross_product( u0, X );
        R = H - L;

        if( !abs_equal( H, L, tol ) )
          return false;
      }
  }
  else
  {
    for( j = 0; j <= nv; j++ )
      for( i = 0; i <= nu; i++ )
      {
        X = euclid( Pw[j][i] );
        R = X - A;

        if( !abs_equal( X, A, tol ) )
          return false;
      }
  }
  return true;
}
// template instantiation
template bool IsPlanar( 
          const int nu, const int nv, const Ptr< Ptr< hpoint<float> > >& Pw,
          const float tol );
template bool IsPlanar( 
          const int nu, const int nv, const Ptr< Ptr< hpoint<double> > >& Pw,
          const double tol );

template bool IsPlanar( 
          const int nu, const int nv, const Ptr< Ptr< point<float> > >& Pw,
          const float tol );
template bool IsPlanar( 
          const int nu, const int nv, const Ptr< Ptr< point<double> > >& Pw,
          const double tol );

/*-------------------------------------------------------------------------
  check if a polyline forms a rectangle
-------------------------------------------------------------------------*/

template< class T >
GULAPI bool isRectangle( int nP, const Ptr< point2<T> >& P,
                         T& u1, T& u2, T& v1, T& v2 )
{
  if( nP != 5 ) return false;

  // last point must be the same as the first one 
  if( (P[0].x != P[4].x) || (P[0].y != P[4].y) )
    return false;

  if( P[0].y == P[1].y )
  {
    if( (P[2].y != P[3].y) || 
        (P[3].x != P[0].x) || (P[2].x != P[1].x) )
      return false;

    if( P[3].y < P[0].y )
    {
      v1 = P[3].y;
      v2 = P[0].y;
    }
    else
    {
      v1 = P[0].y;
      v2 = P[3].y;
    }
    if( P[2].x < P[0].x )
    {
      u1 = P[2].x;
      u2 = P[0].x;
    }
    else
    {
      u1 = P[0].x;
      u2 = P[2].x;
    }
  }
  else if( P[0].x == P[1].x )
  {
    if( (P[2].x != P[3].x) || 
        (P[3].y != P[0].y) || (P[2].y != P[1].y) )
      return false;

    if( P[3].x < P[0].x )
    {
      u1 = P[3].x;
      u2 = P[0].x;
    }
    else
    {
      u1 = P[0].x;
      u2 = P[3].x;
    }
    if( P[2].y < P[0].y )
    {
      u1 = P[2].x;
      u2 = P[0].x;
    }
    else
    {
      u1 = P[0].x;
      u2 = P[2].x;
    }
  }
  else
    return false;

  return true;
}
// template instantiation
template GULAPI bool isRectangle( int nP, const Ptr< point2<float> >& P,
                  float& u1, float& u2, float& v1, float& v2 );
template GULAPI bool isRectangle( int nP, const Ptr< point2<double> >& P,
                  double& u1, double& u2, double& v1, double& v2 );


}

---