lib/cpm/Algorithms/Simplex/DeformationFilter.hxx

   1 #ifndef __CPM__ALGORITHMS__SIMPLEX__DEFORMATIONFILTER__HXX__
   2 #define __CPM__ALGORITHMS__SIMPLEX__DEFORMATIONFILTER__HXX__
   3
   4 // -------------------------------------------------------------------------
   5 template< class M >
   6 void cpm::Algorithms::Simplex::DeformationFilter< M >::
   7 AddForce( TForceFunction* f, const TScalar& w )
   8 {
   9   this->m_Forces.push_back( TFunctionPair( f, w ) );
  10   this->Modified( );
  11 }
  12
  13 // -------------------------------------------------------------------------
  14 template< class M >
  15 cpm::Algorithms::Simplex::DeformationFilter< M >::
  16 DeformationFilter( )
  17   : Superclass( ),
  18     m_Damping( TScalar( 1 ) )
  19 {
  20 }
  21
  22 // -------------------------------------------------------------------------
  23 template< class M >
  24 cpm::Algorithms::Simplex::DeformationFilter< M >::
  25 ~DeformationFilter( )
  26 {
  27 }
  28
  29 // -------------------------------------------------------------------------
  30 template< class M >
  31 void cpm::Algorithms::Simplex::DeformationFilter< M >::
  32 GenerateData( )
  33 {
  34   typedef typename M::VectorType TVector;
  35
  36   const M* in = this->GetInput( );
  37   M* out = this->GetOutput( );
  38   unsigned int nPoints = in->GetNumberOfPoints( );
  39   std::vector< TVector > final_forces[ 2 ];
  40   unsigned int aForces = 0;
  41   final_forces[ 0 ].resize( nPoints, TVector( TScalar( 0 ) ) );
  42   final_forces[ 1 ].resize( nPoints, TVector( TScalar( 0 ) ) );
  43
  44   // Configure forces
  45   typename TFunctionContainer::iterator fIt = this->m_Forces.begin( );
  46   for( ; fIt != this->m_Forces.end( ); ++fIt )
  47     fIt->first->SetMesh( in );
  48
  49   unsigned int leo = 0;
  50   while( leo < 20000 ) // TODO: !!!!
  51   {
  52     leo++;
  53     TScalar damping = TScalar( 1 ) - this->m_Damping;
  54     unsigned int pForces = ( aForces + 1 ) % 2;
  55
  56     // Compute forces
  57     for( TPointId pId = 0; pId < nPoints; ++pId )
  58     {
  59       if( in->FindEdge( pId )->IsOriginInternal( ) )
  60       {
  61         fIt = this->m_Forces.begin( );
  62         TVector f = fIt->first->Evaluate( pId ) * fIt->second;
  63         for( ++fIt; fIt != this->m_Forces.end( ); ++fIt )
  64           f += fIt->first->Evaluate( pId ) * fIt->second;
  65         f += final_forces[ pForces ][ pId ] * damping;
  66         final_forces[ aForces ][ pId ] = f;
  67       }
  68       else
  69         final_forces[ aForces ][ pId ].Fill( TScalar( 0 ) );
  70
  71     } // rof
  72
  73     // Move points
  74     for( TPointId pId = 0; pId < nPoints; ++pId )
  75     {
  76       typename M::PointType pnt = in->GetPoint( pId );
  77       pnt += final_forces[ aForces ][ pId ];
  78       out->SetPoint( pId, pnt );
  79       this->InvokeEvent( PointUpdatedEvent( pId ) );
  80
  81     } // rof
  82     aForces = pForces;
  83
  84     this->InvokeEvent( AllMeshUpdatedEvent( ) );
  85
  86   } //elihw
  87 }
  88
  89 #endif // __CPM__ALGORITHMS__SIMPLEX__DEFORMATIONFILTER__HXX__
  90
  91 // eof - $RCSfile$