Merge branch 'master' of /home/dbv/git/nice

core/algebra/progs/testILSConjugateGradients.cpp

@@ -26,11 +26,23 @@ int main(int argc, char* argv[])
 	std::string logfilename;
 	
 	if ( argc < 2 )
-	  logfilename = "/home/bodesheim/testILS-CGM.log";
-	else
+  {
+    std::cerr << "Warning: you have to specify a log-file with write-access. Program will be closed." << std::endl;
+    return -1;
+  }
+  else
+  {
 	  logfilename = argv[1];
+  }
 	
 	logfile = fopen(logfilename.c_str(), "w");
+  
+  //is the given logfile writable
+  if (logfile == NULL) 
+  {
+    std::cerr << "Error opening file" << std::endl;
+    return -1;
+  }
 	
 	// generate matrix A
 	Matrix A(mySize,mySize,0.0);

core/basics/StringTools.h

@@ -1,4 +1,4 @@
-/** 
+/**
  * @file StringTools.h
  * @brief some string handling routines
  * @author Erik Rodner
@@ -21,82 +21,82 @@ namespace NICE
 class StringTools
 {
 public:
-	/** normalize the std::string value, remove leading and
-	 finishing whitespaces, tabulars, a finishing semicolon
-	 @param s input std::string which will be modified
-	 */
-	static void normalize_string(std::string & s);
-
-	/** split a std::string 
-	 @param s input string
-	 @param splitChar character which was used for seperation
-	 @param list resulting list of sub strings
-	 */
-	static void split(const std::string & s, char splitChar, std::vector<
-			std::string> & list);
-
-	/** split a std::string and convert the strings to double values
-	 @param input std::string value
-	 @param splitChar character which was used as a seperation
-	 @param x resulting std::vector of double values
-	 */
-	static void splitVector(const std::string & s, char splitChar,
-			NICE::Vector & x);
-
-	/** remove a leading and a finishing character
-	 @param value input std::string value
-	 @param trimChar e.g.
-	 */
-	static void trimbounds(std::string & value, char trimChar);
-
-	/** chomp off the end-of-line characters 
-	 @param s input and output string
-	 */
-	static std::string chomp(std::string s);
-
-	/** whatever */
-	static std::string trim(std::string s, const std::string& drop = " ");
-
-	/** trims the filename to be without directories, similar to basename in linux
-	 @param str input string
-	 @param extension specifies if output is with extension or not
-	 */
-	static std::string baseName(const std::string& str, bool extension = true);
-
-	/** perform regular expression matching and return true if the string \c s
-	 *  matches the regular expression in \c regex (GNU regex style!). 
-	 */
-	static bool regexMatch(const std::string & s, const std::string & regex);
-
-	/** perform regular expression matching and return true if the string \c s
-	 *  matches the regular expression in \c regex (GNU regex style!), 
-	 *  sub-pattern matches are stored in submatches and denoted with (...) in
-	 *  the regular expression.
-	 *  @note submatches[0] is the whole match and submatches[i] is the i'th sub-pattern match
-	 */
-	static bool regexMatch(const std::string & s, const std::string & regex,
-			std::vector<std::string> & submatches);
-
-	/** substitute a regular expression match with \c subst */
-	static bool regexSubstitute(std::string & s, const std::string & regex,
-			const std::string & subst);
-
-	/** convert a string to another type using stringstream */
-	template<class F>
-	static F convert(const std::string & s);
-
-	/** convert a string to another type using stringstream and returning true on success and
-	 * false if there is a parsing error */
-	template<class F>
-	static bool convert(const std::string & s, F & i);
-
-	/** convert another type to a string using stringstream */
-	template<class F>
-	static std::string convertToString(const F & s);
-
-	/** convert another type to a string using stringstream */
-	template<class F>
-	static std::string convertToString(const F & s, uint precision);
+    /** normalize the std::string value, remove leading and
+     finishing whitespaces, tabulars, a finishing semicolon
+     @param s input std::string which will be modified
+     */
+    static void normalize_string(std::string & s);
+
+    /** split a std::string
+     @param s input string
+     @param splitChar character which was used for seperation
+     @param list resulting list of sub strings
+     */
+    static void split(const std::string & s, char splitChar, std::vector<
+                      std::string> & list);
+
+    /** split a std::string and convert the strings to double values
+     @param input std::string value
+     @param splitChar character which was used as a seperation
+     @param x resulting std::vector of double values
+     */
+    static void splitVector(const std::string & s, char splitChar,
+                            NICE::Vector & x);
+
+    /** remove a leading and a finishing character
+     @param value input std::string value
+     @param trimChar e.g.
+     */
+    static void trimbounds(std::string & value, char trimChar);
+
+    /** chomp off the end-of-line characters
+     @param s input and output string
+     */
+    static std::string chomp(std::string s);
+
+    /** whatever */
+    static std::string trim(std::string s, const std::string& drop = " ");
+
+    /** trims the filename to be without directories, similar to basename in linux
+     @param str input string
+     @param extension specifies if output is with extension or not
+     */
+    static std::string baseName(const std::string& str, bool extension = true);
+
+    /** perform regular expression matching and return true if the string \c s
+     *  matches the regular expression in \c regex (GNU regex style!).
+     */
+    static bool regexMatch(const std::string & s, const std::string & regex);
+
+    /** perform regular expression matching and return true if the string \c s
+     *  matches the regular expression in \c regex (GNU regex style!),
+     *  sub-pattern matches are stored in submatches and denoted with (...) in
+     *  the regular expression.
+     *  @note submatches[0] is the whole match and submatches[i] is the i'th sub-pattern match
+     */
+    static bool regexMatch(const std::string & s, const std::string & regex,
+                           std::vector<std::string> & submatches);
+
+    /** substitute a regular expression match with \c subst */
+    static bool regexSubstitute(std::string & s, const std::string & regex,
+                                const std::string & subst);
+
+    /** convert a string to another type using stringstream */
+    template<class F>
+    static F convert(const std::string & s);
+
+    /** convert a string to another type using stringstream and returning true on success and
+     * false if there is a parsing error */
+    template<class F>
+    static bool convert(const std::string & s, F & i);
+
+    /** convert another type to a string using stringstream */
+    template<class F>
+    static std::string convertToString(const F & s);
+
+    /** convert another type to a string using stringstream */
+    template<class F>
+    static std::string convertToString(const F & s, uint precision);
 
 };
 
@@ -105,36 +105,36 @@ public:
 template<class F>
 F StringTools::convert(const std::string & s)
 {
-	std::istringstream is(s);
-	F i;
-	is >> i;
-	return i;
+    std::istringstream is(s);
+    F i;
+    is >> i;
+    return i;
 }
 
 template<class F>
 bool StringTools::convert(const std::string & s, F & i)
 {
-	std::istringstream is(s);
-	is >> i;
-	return (!is.fail());
+    std::istringstream is(s);
+    is >> i;
+    return (!is.fail());
 }
 
 template<class F>
 std::string StringTools::convertToString(const F & s)
 {
-	std::ostringstream os;
-	os << s;
-	return os.str();
+    std::ostringstream os;
+    os << s;
+    return os.str();
 }
 
 template<class F>
 std::string StringTools::convertToString(const F & s, uint precision)
 {
-	std::ostringstream os;
-	os.precision(precision);
-	os.setf(std::ios::fixed,std::ios::floatfield);
-	os << s;
-	return os.str();
+    std::ostringstream os;
+    os.precision(precision);
+    os.setf(std::ios::fixed,std::ios::floatfield);
+    os << s;
+    return os.str();
 }
 
 } // namespace

core/image/ColorImageT.h

@@ -1,4 +1,4 @@
-    /*
+/*
  * NICE-Core - efficient algebra and computer vision methods
  *  - libimage - An image library
  * See file License for license information.

core/image/Convert.h

@@ -121,6 +121,18 @@ template<class P> class ImageT;
     */
     template<class P>
     void imageToPseudoColor ( const NICE::ImageT<P> & src, NICE::ColorImage & dst );
+    
+    /**
+    * Convert the Image \c src into a RGB-ColorImage \c dst using pseudo colors.
+    * @param src source gray image
+    * @param dst result image
+    * @param _min desired lower bound (smaller initial values will be set to zero after conversion)
+    * @param _max desired upper bound (larger initial values will be set to 255 after conversion)
+    * @author Alexander Freytag
+    * @date 17-04-2013 (dd-mm-yyyy)
+    */    
+    template<class P>
+    void imageToPseudoColorWithRangeSpecification ( const NICE::ImageT<P> & src, NICE::ColorImage & img, const double & _min, const double & _max );    
 
     /**
     * Convert the matrix \c src into a RGB-ColorImage \c dst using pseudo colors.

core/image/Convert.tcc

@@ -126,27 +126,67 @@ void matrixToPseudoColor ( const NICE::MatrixT<P> & m, NICE::ColorImage & img )
 }
 
 template<class P>
-void imageToPseudoColor ( const NICE::ImageT<P> & m, NICE::ColorImage & img )
+void imageToPseudoColor ( const NICE::ImageT<P> & src, NICE::ColorImage & img )
 {
-  img.resize ( m.width(), m.height() );
+  img.resize ( src.width(), src.height() );
 
   // determine maximum and minimum value in the source image
   // for appropiate scaling
   double max = - std::numeric_limits<double>::max();
   double min = std::numeric_limits<double>::max();
-  for ( size_t x = 0 ; x < ( size_t ) m.width(); x++ )
-    for ( size_t y = 0 ; y < ( size_t ) m.height() ; y++ )
+  for ( size_t x = 0 ; x < ( size_t ) src.width(); x++ )
+    for ( size_t y = 0 ; y < ( size_t ) src.height() ; y++ )
     {
-      double v = m.getPixel ( x, y );
+      double v = src.getPixel ( x, y );
       if ( v > max ) max = v;
       if ( v < min ) min = v;
     }
 
-  for ( size_t y = 0 ; y < ( size_t ) m.height() ; y++ )
-    for ( size_t x = 0 ; x < ( size_t ) m.width(); x++ )
+  for ( size_t y = 0 ; y < ( size_t ) src.height() ; y++ )
+    for ( size_t x = 0 ; x < ( size_t ) src.width(); x++ )
     {
       // scale the grayvalue
-      double val = ( m.getPixel ( x, y ) - min ) / ( max - min );
+      double val = ( src.getPixel ( x, y ) - min ) / ( max - min );
+      double r, g, b;
+      // determine the RGB values
+      convertToPseudoColor ( val, r, g, b );
+      img.setPixel ( x, y, 0, ( int ) ( r*255 ) );
+      img.setPixel ( x, y, 1, ( int ) ( g*255 ) );
+      img.setPixel ( x, y, 2, ( int ) ( b*255 ) );
+    }
+}
+
+template<class P>
+void imageToPseudoColorWithRangeSpecification ( const NICE::ImageT<P> & src, NICE::ColorImage & img, const double & _min, const double & _max )
+{
+  img.resize ( src.width(), src.height() );
+
+  double max ( _max );
+  double min ( _min );
+
+  //consistency check
+  if (max < min )
+  {
+    min = max;
+    max = _min;
+  }
+  
+  //to avoid division by numerical zero
+  double rangeSpecified ( fabs(max - min) );
+  if ( rangeSpecified < 10e-10 )
+  { 
+    max = min + 10e-10;
+    rangeSpecified = fabs(max - min);
+  }
+
+  for ( size_t y = 0 ; y < ( size_t ) src.height() ; y++ )
+    for ( size_t x = 0 ; x < ( size_t ) src.width(); x++ )
+    {
+      // scale the grayvalue
+      double val = ( src.getPixel ( x, y ) - min );
+      val = std::max( val , min ); // check for lower bound
+      val = std::min( val , max ); //check for upper bound
+      val /= rangeSpecified; //appropriate scaling
       double r, g, b;
       // determine the RGB values
       convertToPseudoColor ( val, r, g, b );

core/image/FilterT.h

@@ -96,6 +96,63 @@ class FilterT
     * @return void
     **/
     static void gradientStrength ( const NICE::ImageT<SrcType> &src, NICE::ImageT<DstType> &dst );
+    
+    
+    /**
+    * Blurs Image \c src into the Image \c dst using a mean filter.
+    * @param src  source gray image
+    * @param size mean mask of size (2*\c size + 1) x (2*\c size + 1)
+    * @param dst  Optional buffer to be used as target.<br>
+    *             Create a new Image if \c dst == NULL.<br>
+    *             If \c dst != NULL then size must be equal to \c src 's size!
+    * @throw ImageException will be thrown if \c dst != NULL and the size of \c src and \c dst is not equal.
+    * @return Pointer to Image
+    */
+    ImageT<DstType> * filterMean ( const ImageT<SrcType>& src, const uint& size, ImageT<DstType>* dst = NULL );     
+    
+    /**
+      * Blurs Image \c src into the Image \c dst using a mean filter. This is the
+      * implementation of the filter with a runtime independent of the filter size.
+      * It is significantly faster for large filter sizes. An analysis can be done with the program
+      * testApproxGaussFilter.
+      * The second template parameter gives the type used to store sums of gray-values. The third template
+      * parameter specifies the data type of the temporary image used to handle the seperability of the mean filter.
+      * Boundary handling is done with a cropped filter if use_filtersize_independent_implementation is set to true, otherwise the borders of the image
+      * are ignored.
+      * @param src  source gray image
+      * @param size mean mask of size (2*\c size + 1) x (2*\c size + 1)
+      * @param dst  Optional buffer to be used as target.<br>
+      *             Create a new Image if \c dst == NULL.<br>
+      *             If \c dst != NULL then size must be equal to \c src 's size!
+      * @throw ImageException will be thrown if \c dst != NULL and the size of \c src and \c dst is not equal.
+      * @return Pointer to Image
+      */
+    ImageT<DstType> * filterMeanLargeFS ( const ImageT<SrcType>& src, const uint& size, ImageT<DstType>* dst = NULL );    
+      
+    /**
+    * @brief Blurs Image \c src into the Image \c dst by approximating a gauss filter with a specific standard deviation through multiple mean filters
+    * The size of the mean filters is automatically calculated
+    * The second template parameter gives the type used to store sums of gray-values. The third template parameter,
+    * which is only necessary if use_filtersize_independent_implementation
+    * is set to true, specifies the data type of the temporary image used to handle the seperability of the mean filter.
+    * Boundary handling is done with a cropped filter if use_filtersize_independent_implementation is set to true, otherwise the borders of the image
+    * are ignored.
+      * @param src  source gray image
+      * @param sigma standard deviation of the gauss filter
+      * @param dst  Optional buffer to be used as target.<br>
+      *             Create a new Image if \c dst == NULL.<br>
+      *             If \c dst != NULL then size must be equal to \c src 's size!
+    * @param use_filtersize_independent_implementation  if this flag is set to true, we use the implementation of the mean filter with a runtime
+    *  independent of the filter size. This is beneficial for large filter sizes. However, for small filter sizes (<10) the standard implementation is faster.
+    *  An analysis can be done with the program testApproxGaussFilter.
+      * @throw ImageException will be thrown if \c dst != NULL and the size of \c src and \c dst is not equal
+      *                       or an invalid size \c size is specified.
+      * @return Pointer to Image
+      */
+    NICE::ImageT<DstType> * filterGaussSigmaApproximate ( const NICE::ImageT<SrcType> &src, double sigma, NICE::ImageT<DstType> *dst = NULL,
+        bool use_filtersize_independent_implementation = true );  
+  
+ 
 };
 
 // float specializations using IPP

core/image/FilterT.tcc

@@ -146,4 +146,175 @@ void FilterT<SrcType, CalcType, DstType>::gradientStrength ( const NICE::ImageT<
   return;
 }
 
+template<class SrcType, class CalcType, class DstType>
+ImageT<DstType> * FilterT<SrcType, CalcType, DstType>::filterMean ( const NICE::ImageT<SrcType>& src, const uint& size, ImageT<DstType>* dst )
+{
+  ImageT<DstType>* result = createResultBuffer ( src, dst );
+
+  int isize = size;
+  int is = isize + 1;
+  int msize = ( 2 * isize + 1 ) * ( 2 * isize + 1 );
+  CalcType sum = 0;
+
+  int xend = src.width() - isize;
+
+  const SrcType* pSrc;
+  DstType* pDst;
+  for ( int y = isize; y < src.height() - isize; ++y ) {
+    sum = 0;
+
+    for ( int j = y - isize; j <= y + isize; ++j ) {
+      pSrc = src.getPixelPointerY ( j );
+      for ( uint i = 0; i <= 2*size; ++i, ++pSrc )
+        sum += *pSrc;
+    }
+
+    pSrc = src.getPixelPointerXY ( is, y );
+    pDst = result->getPixelPointerXY ( isize, y );
+    *pDst = sum / msize;
+    ++pDst;
+
+    for ( int x = is; x < xend; ++x, ++pSrc, ++pDst ) {
+      for ( int i = -isize; i <= isize; ++i ) {
+        sum += * ( pSrc + size + i * src.getStepsize() );
+        sum -= * ( pSrc - is + i * src.getStepsize() );
+      }
+      *pDst = sum / msize;
+    }
+  }
+
+  return result;
+}
+
+template<class SrcType, class CalcType, class DstType>
+ImageT<DstType> * FilterT<SrcType, CalcType, DstType>::filterMeanLargeFS ( const ImageT<SrcType>& src, const uint& size, ImageT<DstType>* dst )
+{
+  ImageT<DstType>* result = createResultBuffer ( src, dst );
+  ImageT<CalcType> tmp ( src.width(), src.height() );
+
+  unsigned int isize = size;
+  unsigned int n = 2 * isize + 1;
+  CalcType sum = 0;
+
+  unsigned int xend = src.width();
+  unsigned int yend = src.height();
+
+  // legend of the following comments
+  // current pixel in dst = x
+  // end of filter mask and current pixel in src = e
+  // the examples use a size=2 filter, is=3, n=5
+  for ( int y = 0; y < src.height(); ++y ) {
+    sum = 0;
+
+    const SrcType* pSrc = src.getPixelPointerY ( y );
+    DstType *pDst = tmp.getPixelPointerY ( y );
+
+    // first position: [eoooooooooooooooooo]
+    // last position:  [-eooooooooooooooooo]
+    for ( unsigned int e = 0; e < size; ++e ) {
+      sum += * ( pSrc );
+      pSrc++;
+    }
+    // first position: [x-eoooooooooooooooo]
+    // last position:  [s-x-eoooooooooooooo]
+    for ( unsigned int e = size; e < n; ++e ) {
+      sum += * ( pSrc );
+      pSrc++;
+      *pDst = sum / ( e + 1 );
+      pDst++;
+    }
+    // first position: [os-x-eooooooooooooo]
+    // last position:  [oooooooooooooos-x-e]
+    for ( unsigned int e = n;e < xend; ++e ) {
+      sum -= * ( pSrc - n );
+      sum += * ( pSrc );
+      pSrc++;
+      *pDst = sum / n;
+      pDst++;
+    }
+    // first position: [ooooooooooooooos-x-]
+    // last position:  [oooooooooooooooos-x]
+    for ( unsigned int x = xend - size; x < xend; ++x ) {
+      sum -= * ( pSrc - ( xend - x + 1 ) );
+      *pDst = sum / ( size + xend - x );
+      pDst++;
+    }
+  }
+
+  // now let us filter along the columns
+  long long linestep = src.rowStepsize() / src.bytedepth();
+
+  for ( unsigned int x = 0; x < xend; ++x ) {
+    sum = 0;
+
+    const CalcType* pSrc = tmp.getPixelPointerXY ( x, 0 );
+    SrcType *pDst = result->getPixelPointerXY ( x, 0 );
+
+    for ( unsigned int e = 0; e < size; ++e ) {
+      sum += * ( pSrc );
+      pSrc += linestep;
+    }
+    for ( unsigned int e = size; e < n; ++e ) {
+      sum += * ( pSrc );
+      pSrc += linestep;
+      *pDst = sum / ( e + 1 );
+      pDst += linestep;
+    }
+    for ( unsigned int e = n;e < yend; ++e ) {
+      sum -= * ( pSrc - n * linestep );
+      sum += * ( pSrc );
+      pSrc += linestep;
+      *pDst = sum / n;
+      pDst += linestep;
+    }
+    for ( unsigned int y = yend - size; y < yend; ++y ) {
+      sum -= * ( pSrc - ( yend - y + 1 ) * linestep );
+      *pDst = sum / ( size + yend - y );
+      pDst += linestep;
+    }
+  }
+
+  return result;
+}    
+
+template<class SrcType, class CalcType, class DstType>
+ImageT<DstType> * FilterT<SrcType, CalcType, DstType>::filterGaussSigmaApproximate ( const NICE::ImageT<SrcType> &src, double sigma, NICE::ImageT<DstType> *dst,
+    bool use_filtersize_independent_implementation )
+{
+  // We use the idea of Wells 1986
+  // Efficient Synthesis of Gaussian Filters by Cascaded Uniform Filters
+  // http://ieeexplore.ieee.org/stamp/stamp.jsp?arnumber=04767776
+  // However, the idea is mainly from the lecture script of Prof. Schukat-Talamazzini.
+  //
+  // 4x mean filter of size (2M+1) x (2M+1) is nearly the same as a
+  // gauss filter with sigma = sqrt(4/3(M^2+M))
+  // stddev of a mean filter of size 2M+1 is 1/3 M (M+1)
+  // after applying it 4 times (convolution) the stddev is 4/3 M (M+1)
+  // After some school math stuff (solving quadratic equations),
+  // we can derive the formula used below.
+  // copy-n-paste for gnuplot: plot [0:10] (sqrt(1+3*x*x)-1)/2
+  ImageT<DstType>* result = createResultBuffer ( src, dst );
+  NICE::ImageT<CalcType> tmp ( src.width(), src.height() );
+
+  int M = ( int ) ( ( sqrt ( 1 + 3 * sigma * sigma ) - 1 ) / 2.0 + 0.5 );
+  if ( M < 1 ) {
+    fthrow ( Exception, "Unable to approximate an Gaussian filter of this small scale (sigma=" << sigma << ")" );
+  }
+
+  if ( use_filtersize_independent_implementation )
+  {
+    this->filterMeanLargeFS ( src, M, &tmp );
+    this->filterMeanLargeFS ( tmp, M, result );
+    this->filterMeanLargeFS ( *result, M, &tmp );
+    this->filterMeanLargeFS ( tmp, M, result );
+  } else {
+    this->filterMean ( src, M, &tmp );
+    this->filterMean ( tmp, M, result );
+    this->filterMean ( *result, M, &tmp );
+    this->filterMean ( tmp, M, result );
+  }
+
+  return result;
+}
+
 } //namespace

core/image/MultiChannelImageT.h

@@ -109,7 +109,7 @@ public:
    * @param channel channel
    * @return P mean value of given area
    **/
-  P getIntegralValue(int ulx, int uly, int lrx, int lry, int channel);
+  P getIntegralValue(int ulx, int uly, int lrx, int lry, int channel) const;
 
   /** convert to ice image */
   void convertToGrey( NICE::Image & img, uint channel = 0, bool normalize = true ) const;

core/image/MultiChannelImageT.tcc

@@ -407,9 +407,9 @@ ColorImage MultiChannelImageT<P>::getColor() const
   for ( int y = 0 ; y < ysize; y++ )
     for ( int x = 0 ; x < xsize ; x++, k++ )
     {
-      img.setPixel( x, y, 0, ( int )( data( 0, k ) ) );
-      img.setPixel( x, y, 1, ( int )( data( 1, k ) ) );
-      img.setPixel( x, y, 2, ( int )( data( 2, k ) ) );
+      img.setPixel( x, y, 0, ( int )( data[0][k] ) );
+      img.setPixel( x, y, 1, ( int )( data[1][k] ) );
+      img.setPixel( x, y, 2, ( int )( data[2][k] ) );
     }
 
   //showImage(img);
@@ -448,7 +448,7 @@ void MultiChannelImageT<P>::calcIntegral( uint channel )
 }
 
 template<class P>
-P MultiChannelImageT<P>::getIntegralValue(int ulx, int uly, int lrx, int lry, int channel)
+P MultiChannelImageT<P>::getIntegralValue(int ulx, int uly, int lrx, int lry, int channel) const
 {
   ulx = std::max(ulx-1, -1);
   ulx = std::min(ulx, xsize-1);

core/vector/Distance.tcc

@@ -17,7 +17,7 @@ template<class T>
 T EuclidianDistance<T>::doCalculate(const VectorT<T>& v1, const VectorT<T>& v2) const {
 
     if(v1.size()!=v2.size())
-        _THROW_EVector("Input vectors must have the same size.");
+        _THROW_EVector("Input vectors must have the same size. v1: ");
 
     T dist = T(0);