Commit b31c9df1 authored by Mathias Goldau's avatar Mathias Goldau

[MERGE]

parents d531d348 66aab904
......@@ -88,6 +88,16 @@ void WMEigenSystem::connectors()
void WMEigenSystem::properties()
{
// for selecting the strategy
boost::shared_ptr< WItemSelection > strategies( new WItemSelection() );
strategies->addItem( "LibEigen", "Eigensystem is computed via libEigen and its SelfAdjointEigenSolver" );
strategies->addItem( "Jacobi", "Self implemented Jacobi iterative eigen decomposition" );
m_strategySelector = m_properties->addProperty( "Strategy", "How the eigen system should be computed",
strategies->getSelectorFirst() );
WPropertyHelper::PC_SELECTONLYONE::addTo( m_strategySelector );
WPropertyHelper::PC_NOTEMPTY::addTo( m_strategySelector );
WModule::properties();
}
......@@ -100,6 +110,7 @@ void WMEigenSystem::moduleMain()
{
m_moduleState.setResetable( true, true );
m_moduleState.add( m_tensorIC->getDataChangedCondition() );
m_moduleState.add( m_strategySelector->getCondition() );
ready();
......@@ -243,6 +254,43 @@ WMEigenSystem::TPVODouble::OutTransmitType const WMEigenSystem::eigenFuncDouble(
return computeEigenSystem( m );
}
WMEigenSystem::TPVODouble::OutTransmitType const WMEigenSystem::eigenSolverDouble( TPVODouble::TransmitType const &input )
{
Eigen::Matrix3d m;
m << input[0], input[1], input[2],
input[1], input[3], input[4],
input[2], input[4], input[5];
return applyEigenSolver( m );
}
WMEigenSystem::TPVOFloat::OutTransmitType const WMEigenSystem::eigenSolverFloat( TPVOFloat::TransmitType const &input )
{
Eigen::Matrix3d m;
m << input[0], input[1], input[2],
input[1], input[3], input[4],
input[2], input[4], input[5];
return applyEigenSolver( m );
}
boost::array< double, 12 > WMEigenSystem::applyEigenSolver( const Eigen::Matrix3d &m ) const
{
Eigen::SelfAdjointEigenSolver< Eigen::Matrix3d > es( m );
const Eigen::Matrix3d &evecs = es.eigenvectors();
const Eigen::Vector3d &evals = es.eigenvalues();
boost::array< double, 12 > result = { { evals( 0 ), evecs( 0, 0 ),
evecs( 0, 1 ),
evecs( 0, 2 ),
evals( 1 ), evecs( 1, 0 ),
evecs( 1, 1 ),
evecs( 1, 2 ),
evals( 2 ), evecs( 2, 0 ),
evecs( 2, 1 ),
evecs( 2, 2 ) } }; // NOLINT curly braces
return result;
}
void WMEigenSystem::resetProgress( std::size_t todo, std::string name )
{
if( m_currentProgress )
......@@ -277,15 +325,35 @@ void WMEigenSystem::resetEigenFunction( boost::shared_ptr< WDataSetDTI > tensors
// create a new one
if( tensors->getValueSet()->getDataType() == W_DT_DOUBLE )
{
m_eigenOperationDouble = boost::shared_ptr< TPVODouble >( new TPVODouble( tensors,
boost::bind( &WMEigenSystem::eigenFuncDouble, this, _1 ) ) );
if( m_strategySelector->get().at( 0 )->getName() == "LibEigen" )
{
m_eigenOperationDouble = boost::shared_ptr< TPVODouble >( new TPVODouble( tensors, boost::bind( &WMEigenSystem::eigenSolverDouble, this, _1 ) ) ); // NOLINT line length
}
else if( m_strategySelector->get().at( 0 )->getName() == "Jacobi" )
{
m_eigenOperationDouble = boost::shared_ptr< TPVODouble >( new TPVODouble( tensors, boost::bind( &WMEigenSystem::eigenFuncDouble, this, _1 ) ) ); // NOLINT line length
}
else
{
WAssert( 0, "Bug: Invalid strategy for eigen decomposition selected!" );
}
m_eigenPool = boost::shared_ptr< WThreadedFunctionBase >( new EigenFunctionTypeDouble( W_AUTOMATIC_NB_THREADS, m_eigenOperationDouble ) );
m_moduleState.add( m_eigenPool->getThreadsDoneCondition() );
}
else if( tensors->getValueSet()->getDataType() == W_DT_FLOAT )
{
m_eigenOperationFloat = boost::shared_ptr< TPVOFloat >( new TPVOFloat( tensors,
boost::bind( &WMEigenSystem::eigenFuncFloat, this, _1 ) ) );
if( m_strategySelector->get().at( 0 )->getName() == "LibEigen" )
{
m_eigenOperationFloat = boost::shared_ptr< TPVOFloat >( new TPVOFloat( tensors, boost::bind( &WMEigenSystem::eigenSolverFloat, this, _1 ) ) ); // NOLINT line length
}
else if( m_strategySelector->get().at( 0 )->getName() == "Jacobi" )
{
m_eigenOperationFloat = boost::shared_ptr< TPVOFloat >( new TPVOFloat( tensors, boost::bind( &WMEigenSystem::eigenFuncFloat, this, _1 ) ) ); // NOLINT line length
}
else
{
WAssert( 0, "Bug: Invalid strategy for eigen decomposition selected!" );
}
m_eigenPool = boost::shared_ptr< WThreadedFunctionBase >( new EigenFunctionTypeFloat( W_AUTOMATIC_NB_THREADS, m_eigenOperationFloat ) );
m_moduleState.add( m_eigenPool->getThreadsDoneCondition() );
}
......
......@@ -31,6 +31,8 @@
#include <osg/Geode>
#include "../../ext/Eigen/Eigen"
#include "../../common/math/WTensorFunctions.h"
#include "../../common/WThreadedFunction.h"
#include "../../dataHandler/WThreadedPerVoxelOperation.h"
......@@ -157,7 +159,7 @@ private:
* The function that computes the eigenvectors from the input tensor field.
*
* \param input A subarray of a valueset that consists of the 6 floats that make up the tensor.
* \return The components of the largest eigenvector and the fa value in a 4-double array.
* \return The complete eigen system as double array with 12 components.
*/
TPVOFloat::OutTransmitType const eigenFuncFloat( TPVOFloat::TransmitType const& input );
......@@ -165,10 +167,37 @@ private:
* The function that computes the eigenvectors from the input tensor field.
*
* \param input A subarray of a valueset that consists of the 6 floats that make up the tensor.
* \return The components of the largest eigenvector and the fa value in a 4-double array.
* \return The complete eigen system as double array with 12 components.
*/
TPVODouble::OutTransmitType const eigenFuncDouble( TPVODouble::TransmitType const& input );
/**
* Computes the eigen system for double input parameters via using applyEigenSolver.
*
* \param input A subarray of a valueset that consists of the 6 floats that make up the tensor.
*
* \return The complete eigen system as double array with 12 components.
*/
TPVODouble::OutTransmitType const eigenSolverDouble( TPVODouble::TransmitType const& input );
/**
* Computes the eigen system for double input parameters via using applyEigenSolver.
*
* \param input A subarray of a valueset that consists of the 6 floats that make up the tensor.
*
* \return The complete eigen system as double array with 12 components.
*/
TPVOFloat::OutTransmitType const eigenSolverFloat( TPVOFloat::TransmitType const& input );
/**
* Copies the eigenvalues and eigenvectors from the libEigen output format into the double array.
*
* \param m The symmetric input matrix where only the lower triangular part is considered.
*
* \return The 12 components of the eigen system.
*/
boost::array< double, 12 > applyEigenSolver( const Eigen::Matrix3d& m ) const;
/**
* Is used by every thread to compute the eigensystem for the given tensor.
*
......@@ -211,6 +240,11 @@ private:
* Indicating current work progress.
*/
boost::shared_ptr< WProgress > m_currentProgress;
/**
* List for selecting the strategy.
*/
WPropSelection m_strategySelector;
};
#endif // WMEIGENSYSTEM_H
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment