diff --git a/ChangeLog.md b/ChangeLog.md index 804c2c96e0..a4f26a22b4 100644 --- a/ChangeLog.md +++ b/ChangeLog.md @@ -24,8 +24,7 @@ - Traits class for containers in order to probe their category at compile time. (Jacques-Olivier Lachaud, [#1079](https://github.com/DGtal-team/DGtal/pull/1079)) - -- Generic set operations for arbitrary containers. You may use + - Generic set operations for arbitrary containers. You may use overloaded operators like &, |, -, ^ on arbitrary containers (list, vector, unordered_set, map, etc). (Jacques-Olivier Lachaud, [#1079](https://github.com/DGtal-team/DGtal/pull/1079)) @@ -69,6 +68,14 @@ ## Changes +- *Configuration* + - Types and classes in helper namespaces ```Z2i``` and ```Z3i``` for + ```StdDefs.h``` header (2D and 3D digital geometry with + computations on 32bit integers) are now explicitly instanciated in + the compiled library. This reduces compilation time when such types + are used. (David Coeurjolly, + [#1117](https://github.com/DGtal-team/DGtal/pull/1117)) + - *DEC Package* - DiscreteExteriorCalculus holds both primal and dual sizes of each cell. Subsequent changes have been made to insertSCell. diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index f7428a4fd0..cf4da95c73 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -12,7 +12,7 @@ include(DGtal/arithmetic/ModuleSRC.txt) include(DGtal/kernel/ModuleSRC.txt) include(DGtal/base/ModuleSRC.txt) include(DGtal/io/ModuleSRC.txt) -include(DGtal/math/ModuleSRC.txt) +include(DGtal/helpers/ModuleSRC.txt) ## Board dependency include(Board/ModuleSRC.txt) diff --git a/src/DGtal/base/ModuleSRC.txt b/src/DGtal/base/ModuleSRC.txt index 758117b05f..2276c45b17 100644 --- a/src/DGtal/base/ModuleSRC.txt +++ b/src/DGtal/base/ModuleSRC.txt @@ -3,6 +3,5 @@ SET(DGTAL_SRC ${DGTAL_SRC} DGtal/base/Bits DGtal/base/Clock DGtal/base/Trace - DGtal/base/OrderedAlphabet DGtal/base/Common) diff --git a/src/DGtal/base/OrderedAlphabet.cpp b/src/DGtal/base/OrderedAlphabet.cpp deleted file mode 100644 index c0ee602ac9..0000000000 --- a/src/DGtal/base/OrderedAlphabet.cpp +++ /dev/null @@ -1,518 +0,0 @@ -/** - * 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 . - * - **/ - -/** - * @file OrderedAlphabet.cpp - * @author Jacques-Olivier Lachaud (\c jacques-olivier.lachaud@univ-savoie.fr ) - * Laboratory of Mathematics (CNRS, UMR 5807), University of Savoie, France - * @author Laurent Provot (\c Laurent.Provot@loria.fr ) - * LORIA (CNRS, UMR 7503), Nancy University, France - * - * @date 2010/07/01 - * - * Implementation of methods defined in OrderedAlphabet.h - * - * This file is part of the DGtal library. - */ - -/////////////////////////////////////////////////////////////////////////////// -#include "DGtal/base/OrderedAlphabet.h" -#include "DGtal/arithmetic/ModuloComputer.h" -/////////////////////////////////////////////////////////////////////////////// - -using namespace std; - -/////////////////////////////////////////////////////////////////////////////// -// class OrderedAlphabet -/////////////////////////////////////////////////////////////////////////////// - -/////////////////////////////////////////////////////////////////////////////// -// Standard services - public : - -/** - * Destructor. - */ -DGtal::OrderedAlphabet::~OrderedAlphabet() -{ - if ( myOrder != 0 ) - delete[ ] myOrder; -} - -/** - * @return the current ordered alphabet. - */ -string -DGtal::OrderedAlphabet::orderedAlphabet() const -{ - char *tbl; - tbl = (char *)malloc((myNb + 1)*sizeof(char)); - for ( unsigned int i = 0; i < myNb; ++i ) - { - tbl[ myOrder[ i ] ] = myFirst + i; - } - tbl[ myNb ] = '\0'; - string s( tbl ); - free(tbl); - return s; -} - -/** - * Shift a0 < a1 < ... < an to a1 < ... < an < a0 - */ -void -DGtal::OrderedAlphabet::shiftLeft() -{ - unsigned int* p = myOrder; - unsigned int* q = myOrder + myNb; - while ( p != q ) - { - unsigned int k = *p; - *p = ( k == 0 ) ? myNb - 1 : k - 1; - ++p; - } -} - -/** - * Shift a0 < a1 < ... < an to an < a0 < ... < an-1 - */ -void -DGtal::OrderedAlphabet::shiftRight() -{ - unsigned int* p = myOrder; - unsigned int* q = myOrder + myNb; - while ( p != q ) - { - unsigned int k = *p + 1; - *p = ( k == myNb ) ? 0 : k; - ++p; - } -} - -/** - * Reverse the order a0 < a1 < ... < an to an < ... < a1 < a0 - */ -void -DGtal::OrderedAlphabet::reverse() -{ - unsigned int* p = myOrder; - unsigned int* q = myOrder + myNb; - while ( p != q ) - { - *p = myNb - 1 - *p; - ++p; - } -} - -/** - * Reverse the order a0 < a1 < ... < an to a3 < a2 < a1 < a0 < an < ... - */ -void -DGtal::OrderedAlphabet::reverseAround12() -{ - unsigned int* p = myOrder; - unsigned int* q = myOrder + myNb; - while ( p != q ) - { - *p = ( myNb + 3 - *p ) % myNb; - ++p; - } -} - - - - -/** - * Gives the first lyndon factor of the word [w] starting at - * position [s] in this alphabet. - * - * @param len (returns) the length of the primitive Lyndon factor - * (which starts at position s). - * - * @param nb (returns) the number of times the Lyndon factor appears. - * - * @param w a word - * @param s the starting index in [w]. - * @param e the index after the end in [w] (s mc( modulo ); - index_t i = s; - index_t j = mc.next( s ); - while ( ( j != e ) && ( lessOrEqual( w[ i ], w[ j ] ) ) ) - { - if ( equal( w[ i ], w[ j ] ) ) - mc.increment( i ); - else - i = s; - mc.increment( j ); - } - len = j >= i ? (size_t) ( j - i ) - : (size_t) ( j + modulo - i ); - if (len == 0) - nb = 0; - else - nb = ( (size_t) ( ( j + modulo - s ) % modulo ) ) / len; -} - - -/** - * @param len (returns) the length of the primitive Lyndon factor - * (which starts at position s). - * - * @param nb (returns) the number of times the Lyndon factor appears. - * - * @param w a word which starts with a1 or a2 at position s. - * @param s the starting index in [w]. - * @param e the index after the end in [w] (s mc( modulo ); - ModuloComputer< Integer >::UnsignedInteger i = s; - index_t j = mc.next( s ); - unsigned int p = 1; - unsigned int q = 2; - while ( ( j != e ) && ( lessOrEqual( w[ i ], w[ j ] ) ) ) - { - if ( equal( w[ i ], w[ j ] ) ) - { - if ( j == mc.cast( s + q - 1 ) ) - q += p; - mc.increment( i ); - } - else - { - if ( ( j != mc.cast( s + q - 1 ) ) || ( order ( w[ j ] ) != 2 ) ) - { - len = j; nb = 0; - return false; - } - unsigned int tmp = p; - p = q; - q += q - tmp; - i = s; - } - mc.increment( j ); - } - len = j >= i ? (size_t) ( j - i ) - : (size_t) ( j + modulo - i ); - nb = ( (size_t) ( ( j + modulo - s ) % modulo ) ) / len; - return true; -} - - -/////////////////////////////////////////////////////////////////////////////// -// Interface - public : - -/** - * Writes/Displays the object on an output stream. - * @param out the output stream where the object is written. - */ -void -DGtal::OrderedAlphabet::selfDisplay ( std::ostream & out ) const -{ - out << "[OrderedAlphabet]"; - out << " " << orderedAlphabet() << endl; -} - -/** - * Checks the validity/consistency of the object. - * @return 'true' if the object is valid, 'false' otherwise. - */ -bool -DGtal::OrderedAlphabet::isValid() const -{ - return true; -} - - - -/////////////////////////////////////////////////////////////////////////////// -// ----------------------- MLP services ------------------------------ - -/** - * Extracts the next edge of the minimum length polygon starting from - * position [s] on the word [w]. The alphabet may be modified - * (reversed or shifted). The output alphabet is some a0 < a1 < a2 < ... - * - * @param nb_a1 (returns) the number of letters a1 in the extracted edge (a1 - * in the output alphabet) - * - * @param nb_a2 (returns) the number of letters a2 in the extracted edge (a2 - * in the output alphabet) - * - * @param w the input (cyclic) word (may be modified in the process). - * - * @param s the starting point in the word (updated). - * - * @param cvx (updates) this boolean is flipped only if a change of - * convexity is detected. - * - * @return the number of letters of the extracted edge. - */ -DGtal::OrderedAlphabet::size_t -DGtal::OrderedAlphabet::nextEdge( size_t & nb_a1, - size_t & nb_a2, - std::string & w, - index_t & s, - bool & cvx ) -{ - ModuloComputer< Integer > mc( (const unsigned int)w.size() ); - size_t l; - size_t len; - size_t nb; - bool inC = duvalPPMod( len, nb, w, s, s ); - if ( ! inC ) - // case : change of convexity - { - // JOL : temporary change of letter w[ s ] - char tmp = w[ s ]; - index_t tmp_s = s; - w[ s ] = letter( 2 ); // a3 - this->reverseAround12(); - cvx = ! cvx; - l = nextEdge( nb_a1, nb_a2, w, s, cvx ); - // JOL : former letter is put back in string. - w[ tmp_s ] = tmp; - } - else if ( ( len == 1 ) && ( order( w[ s ] ) == 1 ) ) - // case u=a1 => Quadrant change - { - this->shiftRight(); - s = mc.cast( s + nb ); - nb_a1 = 0; - nb_a2 = nb - 1; - l = nb; - } - else - { // standard (convex) case. - l = len * nb; - char a2 = letter( 2 ); - nb_a1 = len; - nb_a2 = 0; - index_t ss = s; - s = mc.cast( s + l ); - while ( len != 0 ) - { - if ( w[ ss ] == a2 ) ++nb_a2; - mc.increment( ss ); - --len; - } - nb_a1 -= nb_a2; - nb_a1 *= nb; - nb_a2 *= nb; - } - return l; -} - - -/////////////////////////////////////////////////////////////////////////////// -// Internals - private : - -// // -/////////////////////////////////////////////////////////////////////////////// diff --git a/src/DGtal/base/OrderedAlphabet.h b/src/DGtal/base/OrderedAlphabet.h index 5de7f64ba1..8da3602606 100644 --- a/src/DGtal/base/OrderedAlphabet.h +++ b/src/DGtal/base/OrderedAlphabet.h @@ -48,6 +48,7 @@ #include #include "DGtal/base/Common.h" #include "DGtal/kernel/NumberTraits.h" +#include "DGtal/arithmetic/ModuloComputer.h" ////////////////////////////////////////////////////////////////////////////// namespace DGtal diff --git a/src/DGtal/base/OrderedAlphabet.ih b/src/DGtal/base/OrderedAlphabet.ih index dce9780e31..45cceb0457 100644 --- a/src/DGtal/base/OrderedAlphabet.ih +++ b/src/DGtal/base/OrderedAlphabet.ih @@ -158,3 +158,490 @@ DGtal::operator<< ( std::ostream & out, /////////////////////////////////////////////////////////////////////////////// + +/** + * Destructor. + */ +inline +DGtal::OrderedAlphabet::~OrderedAlphabet() +{ + if ( myOrder != 0 ) + delete[ ] myOrder; +} + +/** + * @return the current ordered alphabet. + */ +inline +std::string +DGtal::OrderedAlphabet::orderedAlphabet() const +{ + char *tbl; + tbl = (char *)malloc((myNb + 1)*sizeof(char)); + for ( unsigned int i = 0; i < myNb; ++i ) + { + tbl[ myOrder[ i ] ] = myFirst + i; + } + tbl[ myNb ] = '\0'; + std::string s( tbl ); + free(tbl); + return s; +} + +/** + * Shift a0 < a1 < ... < an to a1 < ... < an < a0 + */ +inline +void +DGtal::OrderedAlphabet::shiftLeft() +{ + unsigned int* p = myOrder; + unsigned int* q = myOrder + myNb; + while ( p != q ) + { + unsigned int k = *p; + *p = ( k == 0 ) ? myNb - 1 : k - 1; + ++p; + } +} + +/** + * Shift a0 < a1 < ... < an to an < a0 < ... < an-1 + */ +inline +void +DGtal::OrderedAlphabet::shiftRight() +{ + unsigned int* p = myOrder; + unsigned int* q = myOrder + myNb; + while ( p != q ) + { + unsigned int k = *p + 1; + *p = ( k == myNb ) ? 0 : k; + ++p; + } +} + +/** + * Reverse the order a0 < a1 < ... < an to an < ... < a1 < a0 + */ +inline +void +DGtal::OrderedAlphabet::reverse() +{ + unsigned int* p = myOrder; + unsigned int* q = myOrder + myNb; + while ( p != q ) + { + *p = myNb - 1 - *p; + ++p; + } +} + +/** + * Reverse the order a0 < a1 < ... < an to a3 < a2 < a1 < a0 < an < ... + */ +inline +void +DGtal::OrderedAlphabet::reverseAround12() +{ + unsigned int* p = myOrder; + unsigned int* q = myOrder + myNb; + while ( p != q ) + { + *p = ( myNb + 3 - *p ) % myNb; + ++p; + } +} + + + + +/** + * Gives the first lyndon factor of the word [w] starting at + * position [s] in this alphabet. + * + * @param len (returns) the length of the primitive Lyndon factor + * (which starts at position s). + * + * @param nb (returns) the number of times the Lyndon factor appears. + * + * @param w a word + * @param s the starting index in [w]. + * @param e the index after the end in [w] (s mc( modulo ); + index_t i = s; + index_t j = mc.next( s ); + while ( ( j != e ) && ( lessOrEqual( w[ i ], w[ j ] ) ) ) + { + if ( equal( w[ i ], w[ j ] ) ) + mc.increment( i ); + else + i = s; + mc.increment( j ); + } + len = j >= i ? (size_t) ( j - i ) + : (size_t) ( j + modulo - i ); + if (len == 0) + nb = 0; + else + nb = ( (size_t) ( ( j + modulo - s ) % modulo ) ) / len; +} + + +/** + * @param len (returns) the length of the primitive Lyndon factor + * (which starts at position s). + * + * @param nb (returns) the number of times the Lyndon factor appears. + * + * @param w a word which starts with a1 or a2 at position s. + * @param s the starting index in [w]. + * @param e the index after the end in [w] (s mc( modulo ); + ModuloComputer< Integer >::UnsignedInteger i = s; + index_t j = mc.next( s ); + unsigned int p = 1; + unsigned int q = 2; + while ( ( j != e ) && ( lessOrEqual( w[ i ], w[ j ] ) ) ) + { + if ( equal( w[ i ], w[ j ] ) ) + { + if ( j == mc.cast( s + q - 1 ) ) + q += p; + mc.increment( i ); + } + else + { + if ( ( j != mc.cast( s + q - 1 ) ) || ( order ( w[ j ] ) != 2 ) ) + { + len = j; nb = 0; + return false; + } + unsigned int tmp = p; + p = q; + q += q - tmp; + i = s; + } + mc.increment( j ); + } + len = j >= i ? (size_t) ( j - i ) + : (size_t) ( j + modulo - i ); + nb = ( (size_t) ( ( j + modulo - s ) % modulo ) ) / len; + return true; +} + + +/////////////////////////////////////////////////////////////////////////////// +// Interface - public : + +/** + * Writes/Displays the object on an output stream. + * @param out the output stream where the object is written. + */ +inline +void +DGtal::OrderedAlphabet::selfDisplay ( std::ostream & out ) const +{ + out << "[OrderedAlphabet]"; + out << " " << orderedAlphabet() << std::endl; +} + +/** + * Checks the validity/consistency of the object. + * @return 'true' if the object is valid, 'false' otherwise. + */ +inline +bool +DGtal::OrderedAlphabet::isValid() const +{ + return true; +} + + + +/////////////////////////////////////////////////////////////////////////////// +// ----------------------- MLP services ------------------------------ + +/** + * Extracts the next edge of the minimum length polygon starting from + * position [s] on the word [w]. The alphabet may be modified + * (reversed or shifted). The output alphabet is some a0 < a1 < a2 < ... + * + * @param nb_a1 (returns) the number of letters a1 in the extracted edge (a1 + * in the output alphabet) + * + * @param nb_a2 (returns) the number of letters a2 in the extracted edge (a2 + * in the output alphabet) + * + * @param w the input (cyclic) word (may be modified in the process). + * + * @param s the starting point in the word (updated). + * + * @param cvx (updates) this boolean is flipped only if a change of + * convexity is detected. + * + * @return the number of letters of the extracted edge. + */ +inline +DGtal::OrderedAlphabet::size_t +DGtal::OrderedAlphabet::nextEdge( size_t & nb_a1, + size_t & nb_a2, + std::string & w, + index_t & s, + bool & cvx ) +{ + ModuloComputer< Integer > mc( (const unsigned int)w.size() ); + size_t l; + size_t len; + size_t nb; + bool inC = duvalPPMod( len, nb, w, s, s ); + if ( ! inC ) + // case : change of convexity + { + // JOL : temporary change of letter w[ s ] + char tmp = w[ s ]; + index_t tmp_s = s; + w[ s ] = letter( 2 ); // a3 + this->reverseAround12(); + cvx = ! cvx; + l = nextEdge( nb_a1, nb_a2, w, s, cvx ); + // JOL : former letter is put back in string. + w[ tmp_s ] = tmp; + } + else if ( ( len == 1 ) && ( order( w[ s ] ) == 1 ) ) + // case u=a1 => Quadrant change + { + this->shiftRight(); + s = mc.cast( s + nb ); + nb_a1 = 0; + nb_a2 = nb - 1; + l = nb; + } + else + { // standard (convex) case. + l = len * nb; + char a2 = letter( 2 ); + nb_a1 = len; + nb_a2 = 0; + index_t ss = s; + s = mc.cast( s + l ); + while ( len != 0 ) + { + if ( w[ ss ] == a2 ) ++nb_a2; + mc.increment( ss ); + --len; + } + nb_a1 -= nb_a2; + nb_a1 *= nb; + nb_a2 *= nb; + } + return l; +} + + +/////////////////////////////////////////////////////////////////////////////// +// Internals - private : + +// // +/////////////////////////////////////////////////////////////////////////////// diff --git a/src/DGtal/geometry/volumes/distance/ExactPredicateLpPowerSeparableMetric.ih b/src/DGtal/geometry/volumes/distance/ExactPredicateLpPowerSeparableMetric.ih index c61352ca75..988002d05e 100644 --- a/src/DGtal/geometry/volumes/distance/ExactPredicateLpPowerSeparableMetric.ih +++ b/src/DGtal/geometry/volumes/distance/ExactPredicateLpPowerSeparableMetric.ih @@ -300,20 +300,6 @@ DGtal::ExactPredicateLpPowerSeparableMetric::closestPower (const Point &o return ClosestBOTH; } //------------------------------------------------------------------------------ -template -inline -typename DGtal::ExactPredicateLpPowerSeparableMetric::Abscissa -DGtal::ExactPredicateLpPowerSeparableMetric::binarySearchHidden(const Abscissa &/*udim*/, - const Abscissa &/*vdim*/, - const Promoted &/*nu*/, - const Promoted &/*nv*/, - const Abscissa &/*lower*/, - const Abscissa &/*upper*/) const -{ - ASSERT(false && "Not Necessary for l_2"); - -} -//------------------------------------------------------------------------------ template inline bool diff --git a/src/DGtal/geometry/volumes/distance/ExactPredicateLpSeparableMetric.ih b/src/DGtal/geometry/volumes/distance/ExactPredicateLpSeparableMetric.ih index 9d970f2bc1..842de45071 100644 --- a/src/DGtal/geometry/volumes/distance/ExactPredicateLpSeparableMetric.ih +++ b/src/DGtal/geometry/volumes/distance/ExactPredicateLpSeparableMetric.ih @@ -318,7 +318,7 @@ DGtal::ExactPredicateLpSeparableMetric::binarySearchHidden(const Abscissa const Abscissa &) const { ASSERT(false && "Not Necessary for l_2"); - + return 0; } //------------------------------------------------------------------------------ template diff --git a/src/DGtal/helpers/ModuleSRC.txt b/src/DGtal/helpers/ModuleSRC.txt new file mode 100644 index 0000000000..d179905412 --- /dev/null +++ b/src/DGtal/helpers/ModuleSRC.txt @@ -0,0 +1,4 @@ + +SET(DGTAL_SRC ${DGTAL_SRC} + DGtal/helpers/StdDefs) + diff --git a/src/DGtal/helpers/StdDefs.cpp b/src/DGtal/helpers/StdDefs.cpp new file mode 100644 index 0000000000..569929dcd2 --- /dev/null +++ b/src/DGtal/helpers/StdDefs.cpp @@ -0,0 +1,67 @@ +/** + * 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 . + * + **/ +/** + * @file StdDefs.cpp + * @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 + * @date 2016/01/19 + * + * This file is part of the DGtal library. + */ + +#include "DGtal/helpers/StdDefs.h" + + +////////////////// +/// Z2i Explicit instanciation +////////////////// +template class DGtal::SpaceND<2, DGtal::Z2i::Integer>; +template class DGtal::PointVector<2, DGtal::Z2i::Integer>; +template class DGtal::KhalimskySpaceND< 2, DGtal::Z2i::Integer > ; +template class DGtal::HyperRectDomain< DGtal::Z2i::Space > ; +template class DGtal::MetricAdjacency< DGtal::Z2i::Space, 1> ; +template class DGtal::MetricAdjacency< DGtal::Z2i::Space, 2> ; +template class DGtal::DigitalTopology< DGtal::Z2i::Adj4, DGtal::Z2i::Adj8 > ; +template class DGtal::DigitalTopology< DGtal::Z2i::Adj8, DGtal::Z2i::Adj4 > ; +template class DGtal::DigitalSetByAssociativeContainer >; +template class DGtal::GridCurve ; +template class DGtal::ExactPredicateLpSeparableMetric ; +template class DGtal::ExactPredicateLpSeparableMetric ; +template class DGtal::ExactPredicateLpPowerSeparableMetric ; +template class DGtal::ExactPredicateLpPowerSeparableMetric ; + + +////////////////// +/// Z3i Explicit instanciation +////////////////// +template class DGtal::SpaceND<3, DGtal::Z3i::Integer>; +template class DGtal::PointVector<3, DGtal::Z3i::Integer>; +template class DGtal::KhalimskySpaceND< 3, DGtal::Z3i::Integer > ; +template class DGtal::HyperRectDomain< DGtal::Z3i::Space > ; +template class DGtal::MetricAdjacency< DGtal::Z3i::Space, 1> ; +template class DGtal::MetricAdjacency< DGtal::Z3i::Space, 2> ; +template class DGtal::MetricAdjacency< DGtal::Z3i::Space, 3> ; +template class DGtal::DigitalTopology< DGtal::Z3i::Adj6, DGtal::Z3i::Adj26 > ; +template class DGtal::DigitalTopology< DGtal::Z3i::Adj26, DGtal::Z3i::Adj6 > ; +template class DGtal::DigitalTopology< DGtal::Z3i::Adj18, DGtal::Z3i::Adj26 > ; +template class DGtal::DigitalTopology< DGtal::Z3i::Adj26, DGtal::Z3i::Adj18 > ; +template class DGtal::DigitalSetByAssociativeContainer >; +template class DGtal::ExactPredicateLpSeparableMetric ; +template class DGtal::ExactPredicateLpSeparableMetric ; +template class DGtal::ExactPredicateLpPowerSeparableMetric ; +template class DGtal::ExactPredicateLpPowerSeparableMetric ; diff --git a/src/DGtal/helpers/StdDefs.h b/src/DGtal/helpers/StdDefs.h index a8f165aa39..b2c50baf03 100644 --- a/src/DGtal/helpers/StdDefs.h +++ b/src/DGtal/helpers/StdDefs.h @@ -183,8 +183,6 @@ namespace DGtal static const DT6_26 dt6_26 = DT6_26( adj6, adj26, JORDAN_DT ); static const DT26_6 dt26_6 = DT26_6( adj26, adj6, JORDAN_DT ); - typedef GridCurve Curve; - typedef ExactPredicateLpSeparableMetric L2Metric; typedef ExactPredicateLpSeparableMetric L1Metric; typedef ExactPredicateLpPowerSeparableMetric L2PowerMetric; diff --git a/src/DGtal/math/AngleLinearMinimizer.cpp b/src/DGtal/math/AngleLinearMinimizer.cpp deleted file mode 100644 index e51e022629..0000000000 --- a/src/DGtal/math/AngleLinearMinimizer.cpp +++ /dev/null @@ -1,522 +0,0 @@ -/** - * 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 . - * - **/ - -/** - * @file AngleLinearMinimizer.cpp - * - * @author Jacques-Olivier Lachaud (\c jacques-olivier.lachaud@univ-savoie.fr ) - * Laboratory of Mathematics (CNRS, UMR 5807), University of Savoie, France - * - * @author (backported from ImaGene by) Bertrand Kerautret (\c kerautre@loria.fr ) - * LORIA (CNRS, UMR 7503), University of Nancy, France - * - * @date 2011/08/31 - * - * Implementation of methods defined in AngleLinearMinimizer.h - * - * This file is part of the DGtal library. - */ - -/////////////////////////////////////////////////////////////////////////////// -#include "DGtal/math/AngleLinearMinimizer.h" -/////////////////////////////////////////////////////////////////////////////// - -/////////////////////////////////////////////////////////////////////////////// -// class AngleLinearMinimizer -/////////////////////////////////////////////////////////////////////////////// - -/////////////////////////////////////////////////////////////////////////////// -// Standard services - public : - - -DGtal::AngleLinearMinimizer::~AngleLinearMinimizer() -{ - reset(); -} - - -DGtal::AngleLinearMinimizer::AngleLinearMinimizer() - : myValues( 0 ) -{ -} - - -void -DGtal::AngleLinearMinimizer::reset() -{ - if ( myValues != 0 ) - { - delete[] myValues; - myValues = 0; - } - mySum = 0.0; - myMax = 0.0; -} - - -void -DGtal::AngleLinearMinimizer::init( unsigned int nbMax ) -{ - reset(); - myValues = new ValueInfo[ nbMax ]; - myMaxSize = nbMax; - mySize = nbMax; - myIsCurveOpen = false; -} - - - - -/////////////////////////////////////////////////////////////////////////////// -// ------------------------- Optimization services -------------------------- - -double -DGtal::AngleLinearMinimizer::getEnergy( unsigned int i1, unsigned int i2 ) const -{ - ModuloComputer mc( size() ); - double E = 0.0; - for ( unsigned int i = mc.next( i1 ); i != i2; ) - { - unsigned int inext = mc.next( i ); - const ValueInfo & vi = this->ro( i ); - const ValueInfo & viprev = this->ro( mc.previous( i ) ); - double dev = AngleComputer::deviation( vi.value, viprev.value ); - E += (dev * dev) / viprev.distToNext; - i = inext; - } - return E; -} - - - -double -DGtal::AngleLinearMinimizer::getFormerEnergy( unsigned int i1, unsigned int i2 ) const -{ - ModuloComputer mc( size() ); - double E = 0.0; - for ( unsigned int i = mc.next( i1 ); i != i2; ) - { - unsigned int inext = mc.next( i ); - const ValueInfo & vi = this->ro( i ); - const ValueInfo & viprev = this->ro( mc.previous( i ) ); - double dev = AngleComputer::deviation( vi.oldValue, viprev.oldValue ); - E += (dev * dev) / viprev.distToNext; - i = inext; - } - return E; -} - - - - -std::vector -DGtal::AngleLinearMinimizer::getGradient() const -{ - ModuloComputer mc( size() ); - - std::vector grad( size() ); - for ( unsigned int i = 0; i < size(); i++ ) - { - unsigned int iprev = mc.previous( i ); - unsigned int inext = mc.next( i ); - const ValueInfo & vi = this->ro( i ); - double val = vi.value; - if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) - { // free extremity to the front/right. - const ValueInfo & viprev = this->ro( iprev ); - double valp = viprev.value; - grad[ i ] = 2.0 * AngleComputer::deviation( val, valp ) / viprev.distToNext; - } - else if ( myIsCurveOpen && ( i == 0 ) ) - { // free extremity to the back/left. - const ValueInfo & vinext = this->ro( inext ); - double valn = vinext.value; - grad[ i ] = -2.0 * AngleComputer::deviation( valn, val ) / vi.distToNext; - } - else - { // standard case. - const ValueInfo & viprev = this->ro( iprev ); - const ValueInfo & vinext = this->ro( inext ); - double valp = viprev.value; - double valn = vinext.value; - grad[ i ] = 2.0*( AngleComputer::deviation( val, valp ) / viprev.distToNext - - AngleComputer::deviation( valn, val ) / vi.distToNext ) ; - } - } - return grad; -} - - - - -std::vector -DGtal::AngleLinearMinimizer::getFormerGradient() const -{ - ModuloComputer< int> mc( size() ); - - std::vector grad( size() ); - for ( unsigned int i = 0; i < size(); i++ ) - { - unsigned int iprev = mc.previous( i ); - unsigned int inext = mc.next( i ); - const ValueInfo & vi = this->ro( i ); - double val = vi.oldValue; - if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) - { // free extremity to the front/right. - const ValueInfo & viprev = this->ro( iprev ); - double valp = viprev.oldValue; - grad[ i ] = 2.0 * AngleComputer::deviation( val, valp ) / viprev.distToNext; - } - else if ( myIsCurveOpen && ( i == 0 ) ) - { // free extremity to the back/left. - const ValueInfo & vinext = this->ro( inext ); - double valn = vinext.oldValue; - grad[ i ] = -2.0 * AngleComputer::deviation( valn, val ) / vi.distToNext; - } - else - { // standard case. - const ValueInfo & viprev = this->ro( iprev ); - const ValueInfo & vinext = this->ro( inext ); - double valp = viprev.oldValue; - double valn = vinext.oldValue; - grad[ i ] = 2.0*( AngleComputer::deviation( val, valp ) / viprev.distToNext - - AngleComputer::deviation( valn, val ) / vi.distToNext ) ; - } - } - return grad; -} - - -double -DGtal::AngleLinearMinimizer::optimize() -{ - return optimize( 0, 0 ); -} - - -double -DGtal::AngleLinearMinimizer::optimize( unsigned int i1, unsigned int i2 ) -{ - ASSERT( size() > 2 ); - ModuloComputer< int> mc( size() ); - - unsigned int i = i1; - // (I) first pass to work on old values. - do - { - ValueInfo & vi = this->rw( i ); - vi.oldValue = vi.value; - // go to next. - i = mc.next( i ); - } - while ( i != i2 ); - this->oneStep( i1, i2 ); - - mySum = 0.0; - myMax = 0.0; - i = i1; - do - { - const ValueInfo & vi = this->ro( i ); - double diff = fabs( AngleComputer::deviation( vi.value, vi.oldValue ) ); - if ( diff > myMax ) - myMax = diff; - mySum += diff; - i = mc.next( i ); - } - while ( i != i2 ); - - return mySum; -} - - -double -DGtal::AngleLinearMinimizer::lastDelta() const -{ - return max(); -} - - -void -DGtal::AngleLinearMinimizer::oneStep( unsigned int i1, unsigned int i2 ) -{ - ModuloComputer< int> mc( size() ); - - double mid; - unsigned int i = i1; - unsigned int iprev = mc.previous( i ); - do - { - unsigned int inext = mc.next( i ); - ValueInfo & vi = this->rw( i ); - if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) - { // free extremity to the front/right. - const ValueInfo & viprev = this->ro( iprev ); - mid = viprev.oldValue; // - viprev.delta; - } - else if ( myIsCurveOpen && ( i == 0 ) ) - { // free extremity to the back/left. - const ValueInfo & vinext = this->ro( inext ); - mid = vinext.oldValue; // + vi.delta; - } - else - { // standard case. - const ValueInfo & viprev = this->ro( iprev ); - const ValueInfo & vinext = this->ro( inext ); - double valp = viprev.oldValue; // - viprev.delta; - double valn = vinext.oldValue; // + vi.delta; - - // old - double y = AngleComputer::deviation( valn, valp ); - mid = ( viprev.distToNext * y ) - / ( vi.distToNext + viprev.distToNext ); - mid = AngleComputer::cast( mid + valp ); - - } - if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; - if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; - double mvt = AngleComputer::deviation( mid, vi.oldValue ); - vi.value = AngleComputer::cast( vi.oldValue + 0.5 * mvt ); - // go to next. - iprev = i; - i = inext; - } - while ( i != i2 ); - - -} - - - -void -DGtal::AngleLinearMinimizerByRelaxation::oneStep( unsigned int i1, unsigned int i2 ) -{ - ModuloComputer< int> mc( size() ); - - double mid; - unsigned int i = i1; - unsigned int iprev = mc.previous( i ); - do - { - unsigned int inext = mc.next( i ); - ValueInfo & vi = this->rw( i ); - // vi.oldValue = vi.value; - if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) - { // free extremity to the front/right. - const ValueInfo & viprev = this->ro( iprev ); - mid = viprev.value; // - viprev.delta; - } - else if ( myIsCurveOpen && ( i == 0 ) ) - { // free extremity to the back/left. - const ValueInfo & vinext = this->ro( inext ); - mid = vinext.oldValue; // + vi.delta; - } - else - { // standard case. - const ValueInfo & viprev = this->ro( iprev ); - const ValueInfo & vinext = this->ro( inext ); - double valp = viprev.value; // - viprev.delta; - double valn = vinext.value; - - // old - double y = AngleComputer::deviation( valn, valp ); - mid = ( viprev.distToNext * y ) - / ( vi.distToNext + viprev.distToNext ); - mid = AngleComputer::cast( mid + valp ); - } - if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; - if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; - vi.value = mid; - iprev = i; - i = inext; - } - while ( i != i2 ); -} - - - - -double -DGtal::AngleLinearMinimizerByRelaxation::lastDelta() const -{ - return max(); -} - - - -void -DGtal::AngleLinearMinimizerByGradientDescent::oneStep( unsigned int i1, unsigned int i2 ) -{ - ModuloComputer< int> mc( size() ); - - std::vector grad ( getFormerGradient() ); - double mid; - unsigned int i = i1; - do - { - unsigned int inext = mc.next( i ); - ValueInfo & vi = this->rw( i ); - double val = vi.oldValue; - mid = AngleComputer::cast( val - myStep*grad[ i ] ); - if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; - if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; - vi.value = mid; - // go to next. - i = inext; - } - while ( i != i2 ); -} - - - -double -DGtal::AngleLinearMinimizerByGradientDescent::lastDelta() const -{ - std::vector grad ( getFormerGradient() ); - double ninf = 0.0; - for ( unsigned int i = 0; i < size(); i++ ) - { - const ValueInfo & vi = this->ro( i ); - if ( vi.value != vi.oldValue ) - { - double n = fabs( grad[ i ] ); - if ( n > ninf ) ninf = n; - } - } - return ninf; -} - - - - -void -DGtal::AngleLinearMinimizerByAdaptiveStepGradientDescent::oneStep( unsigned int i1, unsigned int i2 ) -{ - ModuloComputer< int> mc( size() ); - std::vector grad ( getFormerGradient() ); - - double mid; - unsigned int i = i1; - do - { - unsigned int inext = mc.next( i ); - ValueInfo & vi = this->rw( i ); - double val = vi.oldValue; - mid = AngleComputer::cast( val - myStep*grad[ i ] ); - if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; - if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; - vi.value = mid; - // go to next. - i = inext; - } - while ( i != i2 ); - - double E1 = getFormerEnergy( i1, i2 ); - double E2 = getEnergy( i1, i2 ); - if ( E1 <= E2 ) - { - myStep /= 4.0; - } - else - { - /* doubling step. */ - myStep *= 2.0; - } -} - - - -double -DGtal::AngleLinearMinimizerByAdaptiveStepGradientDescent::lastDelta() const -{ - std::vector grad ( getFormerGradient() ); - double ninf = 0.0; - for ( unsigned int i = 0; i < size(); i++ ) - { - const ValueInfo & vi = this->ro( i ); - if ( vi.value != vi.oldValue ) - { - double n = fabs( grad[ i ] ); - if ( n > ninf ) ninf = n; - } - } - return ninf; -} - - - - -/////////////////////////////////////////////////////////////////////////////// -// Interface - public : - -/** - * Writes/Displays the object on an output stream. - * @param out the output stream where the object is written. - */ -void -DGtal::AngleLinearMinimizer::selfDisplay ( std::ostream & aStream ) const -{ - aStream << "[AngleLinearMinimizer::standard 0.5]"; -} - - -/** - * Writes/Displays the object on an output stream. - * @param aStream the output stream where the object is written. - */ -void -DGtal::AngleLinearMinimizerByRelaxation::selfDisplay( std::ostream & aStream ) const -{ - aStream << "[LinearMinimizer::relaxation]"; -} - -/** - * Writes/Displays the object on an output stream. - * @param aStream the output stream where the object is written. - */ -void -DGtal::AngleLinearMinimizerByGradientDescent::selfDisplay( std::ostream& aStream ) const -{ - aStream << "[LinearMinimizer::gradient descent " << myStep << "]"; -} - -/** - * Writes/Displays the object on an output stream. - * @param aStream the output stream where the object is written. - */ -void -DGtal::AngleLinearMinimizerByAdaptiveStepGradientDescent::selfDisplay( std::ostream& aStream ) const -{ - aStream << "[LinearMinimizer::gradient descent with adaptive step " << myStep << "]"; -} - -/** - * Checks the validity/consistency of the object. - * @return 'true' if the object is valid, 'false' otherwise. - */ -bool -DGtal::AngleLinearMinimizer::isValid() const -{ - return true; -} - - - -/////////////////////////////////////////////////////////////////////////////// -// Internals - private : - -// // -/////////////////////////////////////////////////////////////////////////////// diff --git a/src/DGtal/math/AngleLinearMinimizer.ih b/src/DGtal/math/AngleLinearMinimizer.ih index cd606036db..e15d1246f7 100644 --- a/src/DGtal/math/AngleLinearMinimizer.ih +++ b/src/DGtal/math/AngleLinearMinimizer.ih @@ -192,6 +192,511 @@ DGtal::AngleLinearMinimizer::setIsCurveOpen( bool isCurveOpen) +inline +DGtal::AngleLinearMinimizer::~AngleLinearMinimizer() +{ + reset(); +} + + +inline +DGtal::AngleLinearMinimizer::AngleLinearMinimizer() + : myValues( 0 ) +{ +} + + +inline +void +DGtal::AngleLinearMinimizer::reset() +{ + if ( myValues != 0 ) + { + delete[] myValues; + myValues = 0; + } + mySum = 0.0; + myMax = 0.0; +} + + +inline +void +DGtal::AngleLinearMinimizer::init( unsigned int nbMax ) +{ + reset(); + myValues = new ValueInfo[ nbMax ]; + myMaxSize = nbMax; + mySize = nbMax; + myIsCurveOpen = false; +} + + + + +/////////////////////////////////////////////////////////////////////////////// +// ------------------------- Optimization services -------------------------- + +inline +double +DGtal::AngleLinearMinimizer::getEnergy( unsigned int i1, unsigned int i2 ) const +{ + ModuloComputer mc( size() ); + double E = 0.0; + for ( unsigned int i = mc.next( i1 ); i != i2; ) + { + unsigned int inext = mc.next( i ); + const ValueInfo & vi = this->ro( i ); + const ValueInfo & viprev = this->ro( mc.previous( i ) ); + double dev = AngleComputer::deviation( vi.value, viprev.value ); + E += (dev * dev) / viprev.distToNext; + i = inext; + } + return E; +} + + + +inline +double +DGtal::AngleLinearMinimizer::getFormerEnergy( unsigned int i1, unsigned int i2 ) const +{ + ModuloComputer mc( size() ); + double E = 0.0; + for ( unsigned int i = mc.next( i1 ); i != i2; ) + { + unsigned int inext = mc.next( i ); + const ValueInfo & vi = this->ro( i ); + const ValueInfo & viprev = this->ro( mc.previous( i ) ); + double dev = AngleComputer::deviation( vi.oldValue, viprev.oldValue ); + E += (dev * dev) / viprev.distToNext; + i = inext; + } + return E; +} + + + + +inline +std::vector +DGtal::AngleLinearMinimizer::getGradient() const +{ + ModuloComputer mc( size() ); + + std::vector grad( size() ); + for ( unsigned int i = 0; i < size(); i++ ) + { + unsigned int iprev = mc.previous( i ); + unsigned int inext = mc.next( i ); + const ValueInfo & vi = this->ro( i ); + double val = vi.value; + if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) + { // free extremity to the front/right. + const ValueInfo & viprev = this->ro( iprev ); + double valp = viprev.value; + grad[ i ] = 2.0 * AngleComputer::deviation( val, valp ) / viprev.distToNext; + } + else if ( myIsCurveOpen && ( i == 0 ) ) + { // free extremity to the back/left. + const ValueInfo & vinext = this->ro( inext ); + double valn = vinext.value; + grad[ i ] = -2.0 * AngleComputer::deviation( valn, val ) / vi.distToNext; + } + else + { // standard case. + const ValueInfo & viprev = this->ro( iprev ); + const ValueInfo & vinext = this->ro( inext ); + double valp = viprev.value; + double valn = vinext.value; + grad[ i ] = 2.0*( AngleComputer::deviation( val, valp ) / viprev.distToNext + - AngleComputer::deviation( valn, val ) / vi.distToNext ) ; + } + } + return grad; +} + + + + +inline +std::vector +DGtal::AngleLinearMinimizer::getFormerGradient() const +{ + ModuloComputer< int> mc( size() ); + + std::vector grad( size() ); + for ( unsigned int i = 0; i < size(); i++ ) + { + unsigned int iprev = mc.previous( i ); + unsigned int inext = mc.next( i ); + const ValueInfo & vi = this->ro( i ); + double val = vi.oldValue; + if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) + { // free extremity to the front/right. + const ValueInfo & viprev = this->ro( iprev ); + double valp = viprev.oldValue; + grad[ i ] = 2.0 * AngleComputer::deviation( val, valp ) / viprev.distToNext; + } + else if ( myIsCurveOpen && ( i == 0 ) ) + { // free extremity to the back/left. + const ValueInfo & vinext = this->ro( inext ); + double valn = vinext.oldValue; + grad[ i ] = -2.0 * AngleComputer::deviation( valn, val ) / vi.distToNext; + } + else + { // standard case. + const ValueInfo & viprev = this->ro( iprev ); + const ValueInfo & vinext = this->ro( inext ); + double valp = viprev.oldValue; + double valn = vinext.oldValue; + grad[ i ] = 2.0*( AngleComputer::deviation( val, valp ) / viprev.distToNext + - AngleComputer::deviation( valn, val ) / vi.distToNext ) ; + } + } + return grad; +} + + +inline +double +DGtal::AngleLinearMinimizer::optimize() +{ + return optimize( 0, 0 ); +} + + +inline +double +DGtal::AngleLinearMinimizer::optimize( unsigned int i1, unsigned int i2 ) +{ + ASSERT( size() > 2 ); + ModuloComputer< int> mc( size() ); + + unsigned int i = i1; + // (I) first pass to work on old values. + do + { + ValueInfo & vi = this->rw( i ); + vi.oldValue = vi.value; + // go to next. + i = mc.next( i ); + } + while ( i != i2 ); + this->oneStep( i1, i2 ); + + mySum = 0.0; + myMax = 0.0; + i = i1; + do + { + const ValueInfo & vi = this->ro( i ); + double diff = fabs( AngleComputer::deviation( vi.value, vi.oldValue ) ); + if ( diff > myMax ) + myMax = diff; + mySum += diff; + i = mc.next( i ); + } + while ( i != i2 ); + + return mySum; +} + + +inline +double +DGtal::AngleLinearMinimizer::lastDelta() const +{ + return max(); +} + + +inline +void +DGtal::AngleLinearMinimizer::oneStep( unsigned int i1, unsigned int i2 ) +{ + ModuloComputer< int> mc( size() ); + + double mid; + unsigned int i = i1; + unsigned int iprev = mc.previous( i ); + do + { + unsigned int inext = mc.next( i ); + ValueInfo & vi = this->rw( i ); + if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) + { // free extremity to the front/right. + const ValueInfo & viprev = this->ro( iprev ); + mid = viprev.oldValue; // - viprev.delta; + } + else if ( myIsCurveOpen && ( i == 0 ) ) + { // free extremity to the back/left. + const ValueInfo & vinext = this->ro( inext ); + mid = vinext.oldValue; // + vi.delta; + } + else + { // standard case. + const ValueInfo & viprev = this->ro( iprev ); + const ValueInfo & vinext = this->ro( inext ); + double valp = viprev.oldValue; // - viprev.delta; + double valn = vinext.oldValue; // + vi.delta; + + // old + double y = AngleComputer::deviation( valn, valp ); + mid = ( viprev.distToNext * y ) + / ( vi.distToNext + viprev.distToNext ); + mid = AngleComputer::cast( mid + valp ); + + } + if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; + if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; + double mvt = AngleComputer::deviation( mid, vi.oldValue ); + vi.value = AngleComputer::cast( vi.oldValue + 0.5 * mvt ); + // go to next. + iprev = i; + i = inext; + } + while ( i != i2 ); + + +} + + + +inline +void +DGtal::AngleLinearMinimizerByRelaxation::oneStep( unsigned int i1, unsigned int i2 ) +{ + ModuloComputer< int> mc( size() ); + + double mid; + unsigned int i = i1; + unsigned int iprev = mc.previous( i ); + do + { + unsigned int inext = mc.next( i ); + ValueInfo & vi = this->rw( i ); + // vi.oldValue = vi.value; + if ( myIsCurveOpen && ( i == ( size() - 1 ) ) ) + { // free extremity to the front/right. + const ValueInfo & viprev = this->ro( iprev ); + mid = viprev.value; // - viprev.delta; + } + else if ( myIsCurveOpen && ( i == 0 ) ) + { // free extremity to the back/left. + const ValueInfo & vinext = this->ro( inext ); + mid = vinext.oldValue; // + vi.delta; + } + else + { // standard case. + const ValueInfo & viprev = this->ro( iprev ); + const ValueInfo & vinext = this->ro( inext ); + double valp = viprev.value; // - viprev.delta; + double valn = vinext.value; + + // old + double y = AngleComputer::deviation( valn, valp ); + mid = ( viprev.distToNext * y ) + / ( vi.distToNext + viprev.distToNext ); + mid = AngleComputer::cast( mid + valp ); + } + if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; + if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; + vi.value = mid; + iprev = i; + i = inext; + } + while ( i != i2 ); +} + + + + +inline +double +DGtal::AngleLinearMinimizerByRelaxation::lastDelta() const +{ + return max(); +} + + + +inline +void +DGtal::AngleLinearMinimizerByGradientDescent::oneStep( unsigned int i1, unsigned int i2 ) +{ + ModuloComputer< int> mc( size() ); + + std::vector grad ( getFormerGradient() ); + double mid; + unsigned int i = i1; + do + { + unsigned int inext = mc.next( i ); + ValueInfo & vi = this->rw( i ); + double val = vi.oldValue; + mid = AngleComputer::cast( val - myStep*grad[ i ] ); + if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; + if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; + vi.value = mid; + // go to next. + i = inext; + } + while ( i != i2 ); +} + + + +inline +double +DGtal::AngleLinearMinimizerByGradientDescent::lastDelta() const +{ + std::vector grad ( getFormerGradient() ); + double ninf = 0.0; + for ( unsigned int i = 0; i < size(); i++ ) + { + const ValueInfo & vi = this->ro( i ); + if ( vi.value != vi.oldValue ) + { + double n = fabs( grad[ i ] ); + if ( n > ninf ) ninf = n; + } + } + return ninf; +} + + + + +inline +void +DGtal::AngleLinearMinimizerByAdaptiveStepGradientDescent::oneStep( unsigned int i1, unsigned int i2 ) +{ + ModuloComputer< int> mc( size() ); + std::vector grad ( getFormerGradient() ); + + double mid; + unsigned int i = i1; + do + { + unsigned int inext = mc.next( i ); + ValueInfo & vi = this->rw( i ); + double val = vi.oldValue; + mid = AngleComputer::cast( val - myStep*grad[ i ] ); + if ( AngleComputer::less( mid, vi.min ) ) mid = vi.min; + if ( AngleComputer::less( vi.max, mid ) ) mid = vi.max; + vi.value = mid; + // go to next. + i = inext; + } + while ( i != i2 ); + + double E1 = getFormerEnergy( i1, i2 ); + double E2 = getEnergy( i1, i2 ); + if ( E1 <= E2 ) + { + myStep /= 4.0; + } + else + { + /* doubling step. */ + myStep *= 2.0; + } +} + + + +inline +double +DGtal::AngleLinearMinimizerByAdaptiveStepGradientDescent::lastDelta() const +{ + std::vector grad ( getFormerGradient() ); + double ninf = 0.0; + for ( unsigned int i = 0; i < size(); i++ ) + { + const ValueInfo & vi = this->ro( i ); + if ( vi.value != vi.oldValue ) + { + double n = fabs( grad[ i ] ); + if ( n > ninf ) ninf = n; + } + } + return ninf; +} + + + + +/////////////////////////////////////////////////////////////////////////////// +// Interface - public : + +/** + * Writes/Displays the object on an output stream. + * @param out the output stream where the object is written. + */ +inline +void +DGtal::AngleLinearMinimizer::selfDisplay ( std::ostream & aStream ) const +{ + aStream << "[AngleLinearMinimizer::standard 0.5]"; +} + + +/** + * Writes/Displays the object on an output stream. + * @param aStream the output stream where the object is written. + */ +inline +void +DGtal::AngleLinearMinimizerByRelaxation::selfDisplay( std::ostream & aStream ) const +{ + aStream << "[LinearMinimizer::relaxation]"; +} + +/** + * Writes/Displays the object on an output stream. + * @param aStream the output stream where the object is written. + */ +inline +void +DGtal::AngleLinearMinimizerByGradientDescent::selfDisplay( std::ostream& aStream ) const +{ + aStream << "[LinearMinimizer::gradient descent " << myStep << "]"; +} + +/** + * Writes/Displays the object on an output stream. + * @param aStream the output stream where the object is written. + */ +inline +void +DGtal::AngleLinearMinimizerByAdaptiveStepGradientDescent::selfDisplay( std::ostream& aStream ) const +{ + aStream << "[LinearMinimizer::gradient descent with adaptive step " << myStep << "]"; +} + +/** + * Checks the validity/consistency of the object. + * @return 'true' if the object is valid, 'false' otherwise. + */ +inline +bool +DGtal::AngleLinearMinimizer::isValid() const +{ + return true; +} + + + +/////////////////////////////////////////////////////////////////////////////// +// Internals - private : + +// // +/////////////////////////////////////////////////////////////////////////////// + + + + /////////////////////////////////////////////////////////////////////////////// // Implementation of inline functions and external operators // diff --git a/src/DGtal/math/ModuleSRC.txt b/src/DGtal/math/ModuleSRC.txt deleted file mode 100644 index cf73f5b4e1..0000000000 --- a/src/DGtal/math/ModuleSRC.txt +++ /dev/null @@ -1,6 +0,0 @@ -## Sources associated to the module math -## - -SET(DGTAL_SRC ${DGTAL_SRC} - DGtal/math/AngleLinearMinimizer) -