approximate_min_ellipsoid_d_debug.h

来自「很多二维 三维几何计算算法 C++ 类库」· C头文件 代码 · 共 728 行 · 第 1/2 页

// Copyright (c) 1997-2001  ETH Zurich (Switzerland).
// All rights reserved.
//
// This file is part of CGAL (www.cgal.org); you may redistribute it under
// the terms of the Q Public License version 1.0.
// See the file LICENSE.QPL distributed with CGAL.
//
// Licensees holding a valid commercial license may use this file in
// accordance with the commercial license agreement provided with the software.
//
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
//
// $URL: svn+ssh://scm.gforge.inria.fr/svn/cgal/branches/CGAL-3.3-branch/Approximate_min_ellipsoid_d/include/CGAL/Approximate_min_ellipsoid_d/Approximate_min_ellipsoid_d_debug.h $
// $Id: Approximate_min_ellipsoid_d_debug.h 37176 2007-03-17 08:35:57Z afabri $
// 
//
// Author(s)     : Kaspar Fischer <[e-mail]>

#ifndef CGAL_APPROX_MIN_ELL_D_DEBUG_H
#define CGAL_APPROX_MIN_ELL_D_DEBUG_H

#include <cmath>
#include <map>
#include <iostream>
#include <fstream>

namespace CGAL {

  // We hide the following debugging utilities in a "private" namespace:
  namespace Approximate_min_ellipsoid_d_impl {

    // For debugging only:
    template<typename Iterator>
    void print_matrix(int d, const char *name, Iterator A)
    {
      std::cout << name << ":= Matrix([\n";
      for (int i=0; i<d; ++i) {
	std::cout << "  [ ";
	for (int j=0; j<d; ++j) {
	  std::cout << std::setprecision(30) << A[i+j*d];
	  if (j<d-1)
	    std::cout << ", ";
	}
	std::cout << " ]";
	if (i<d-1)
	  std::cout << ",";
	std::cout << "\n";
      }
      std::cout << "]);\n";	
    }

    // For debugging only:
    template<typename Iterator>
    void print_vector(int d, const char *name, Iterator v)
    {
      std::cout << name << ":= Vector([\n";
      for (int j=0; j<d; ++j) {
	std::cout << std::setprecision(30) << v[j];
	if (j<d-1)
	  std::cout << ", ";
      }
      std::cout << "]);\n";	
    }

    #ifdef CGAL_APPEL_LOG_MODE
    class Logger
    // A singleton class which sends debugging comments to log files.
    // (A class is called a "singleton" if at most one instance of it
    // may ever exist.)  By calling void log(channel,msg) you can send
    // the string msg to the file with name channel.
    {
    private: // private members:
      typedef std::map<std::string,std::ofstream*> Streams;
      Streams channels;           // a collection of pairs (k,v) where
                                  // k is the file-name and v is the
                                  // (open) stream associated with k
  
    private: // (Construction and destruction are private to prevent
  	   // more than one instantiation.)
  
      Logger() {}
      Logger(const Logger&);
      Logger& operator=(const Logger&);
  
      ~Logger()
      {
        // we walk through the list of all opened files and close them:
        for (Streams::iterator it = channels.begin();
  	   it != channels.end(); ++it) {
  	(*it).second->close();
  	delete (*it).second;
        }
      }
      
    public: // access and routines:
  
      static Logger& instance()
      // Returns a reference to the only existing instance of this class:
      {
        // Here's where we maintain the only instance: (Notice that it
        // gets constructed automatically the first time instance() is
        // called, and that it gets disposed of (if ever contructed) at
        // program termination.)
        static Logger instance;
        return instance;
      }
  
      void log(const char *channel,const std::string& msg)
      // If this is the first call to log with string channel as the
      // first parameter, then the file with name channel.log is
      // opened (at the beginning) for writing and msg is written to
      // it.  Otherwise, the string msg is opened to the already open
      // file channel.log.
      {
        const std::string name(channel);
        Streams::iterator it = channels.find(name);
        
        // have we already opened this file?
        if (it != channels.end()) {
  	// If so, then just append the message:
  	*(*it).second << msg;
  	(*it).second->flush();
  	
        } else {
  	// If we haven't seen 'name' before, we create a new file:
  	using std::ofstream;
  	ofstream *o = new ofstream((name+".log").c_str(),
  				   ofstream::out|ofstream::trunc);
  	channels[name] = o;
  	*o << msg;
        }
      }
    };
    #endif // CGAL_APPEL_LOG_MODE

  } // end of namespace Approximate_min_ellipsoid_d_impl
} // end of namespace CGAL

#ifdef CGAL_APPEL_TIMER_MODE
#include <sys/time.h>
#include <sys/resource.h>

namespace CGAL {
  namespace Approximate_min_ellipsoid_d_impl {

    // The following routine is taken from file mptimeval.h from
    // "Matpack Library Release 1.7.1" which is copyright (C) 1991-2002
    // by Berndt M. Gammel.  It works on the timeval struct defined in
    // sys/time.h:
    //
    //       struct timeval {
    //         long tv_sec;        /* seconds */
    //         long tv_usec;       /* microseconds */
    //       };
    //
    inline timeval& operator-=(timeval &t1,const timeval &t2)
    {
      t1.tv_sec -= t2.tv_sec;
      if ( (t1.tv_usec -= t2.tv_usec) < 0 ) {
        --t1.tv_sec;
        t1.tv_usec += 1000000;
      }
      return t1;
    }

    // (Adapted from the above.)
    inline timeval& operator+=(timeval &t1,const timeval &t2)
    {
      t1.tv_sec += t2.tv_sec;
      if ( (t1.tv_usec += t2.tv_usec) > 1000000) {
        ++t1.tv_sec;
        t1.tv_usec -= 1000000;
      }
      return t1;
    }
    
    class Timer
      // A singleton class which maintains a collection of named timers.
      // (A class is called a "singleton" if at most one instance of it
      // may ever exist.)  The following routines are provided:
      //
      // - start(name): If this is the first time start() has been
      //   called with name as the first parameter, then a new timer is
      //   created and started.  Otherwise, the timer with name name is
      //   restarted.
      //
      // - lapse(name): Retuns the number of seconds which have elapsed
      //   since start(name) was called last.
      //   Precondition: start(name) has been called once.
    {
    private: // (Construction and destruction are private to prevent
  	   // more than one instantiation.)
  
      Timer() {}
      Timer(const Timer&);
      Timer& operator=(const Timer&);
  
    public: // access and routines:
  
      static Timer& instance()
      // Returns a reference to the only existing instance of this class:
      {
        // Here's where we maintain the only instance: (Notice that it
        // gets constructed automatically the first time instance() is
        // called, and that it gets disposed of (if ever contructed) at
        // program termination.)
        static Timer instance;
        return instance;
      }
  
      void start(const char *timer_name)
      {
        // fetch current usage:
        rusage now;
        int status = getrusage(RUSAGE_SELF,&now);
        assert(status == 0);
        
        // save it:
        timers[std::string(timer_name)] = now.ru_utime;
      }
  
      float lapse(const char *name)
      {
        // assert that start(name) has been called before:
        assert(timers.find(std::string(name)) != timers.end());
        
        // get current usage:
        rusage now;
        int status = getrusage(RUSAGE_SELF,&now);
        assert(status == 0);
        
        // compute elapsed usage:
        now.ru_utime -= (*timers.find(std::string(name))).second;
        return now.ru_utime.tv_sec + now.ru_utime.tv_usec * 1e-6;
      }
      
    private: // private members:
      typedef std::map<std::string,timeval> Timers;
      Timers timers;              // a collection of pairs (k,v) where
                                  // k is the timer name and v is the
                                  // (started) timer associated with k
    };
  } // end of namespace Approximate_min_ellipsoid_d_impl
} // end of namespace CGAL
#endif // CGAL_APPEL_TIMER_MODE

namespace CGAL {
  namespace Approximate_min_ellipsoid_d_impl {

    template <typename T>
    std::string tostr(const T& t) {
      std::stringstream strm;
      strm << t;
      return strm.str();
    }
    
    class Eps_export_2 {
    // An instance of the following class accepts circles and ellipses
    // and procudes an Enhanced-PostScript figure.
    public:
      enum Stroke_mode { Solid=0, Solid_filled=1, Dashed=2 };
      enum Label_mode { None, Angle, Random_angle };
      
    private:
      std::vector<double> bb;            // bounding box
      Stroke_mode bm;                    // current stroke mode
      Label_mode lm;                     // current label mode
      double next_angle;                 // next angle to use (for Angle mode)
      int count;                         // counts the objects (starting at 1)
      std::string next_label;            // next label to use
      std::ostringstream body;           // buffered output
      bool adjust_bb;
      double zoom;
      std::string filename;
      
    public: // construction and destruction:
      
      // Begins exporting to file filename.  Sets the current stroke mode
      // to Solid and the current label-mode to Random_angle.
      Eps_export_2(std::string filename, double zoom, int /* seed */ = 0) 
	: bb(4,0.0), bm(Solid), lm(Random_angle), count(1), 
	  next_label(default_label(1)), adjust_bb(true),
	  zoom(zoom), filename(filename)
      {}
      
      // Ends exporting and closes the file.
      ~Eps_export_2()
      {
	// open output file:
	using std::endl;
	std::ofstream f(filename.c_str());
	if (!f)
	  std::cerr << "Couldn't open file " << filename << "." << endl;
	
	// write header:
	const double margin = 10.0;
	f << "%!PS-Adobe-2.0 EPSF-2.0" << endl
	  << "%%Title: Some balls" << endl
	  << "%%Creator: a program by kf." << endl
	  << "%%CreationDate: now" << endl
	  << "%%For: you, the unknown" << endl
	  << "%%Pages: 1" << endl
	  << "%%DocumentFonts: Helvetica" << endl
	  << "%%BoundingBox: "
	  << bb[0]*zoom-margin << " "
	  << bb[1]*zoom-margin << " "
	  << bb[2]*zoom+margin << " "
	  << bb[3]*zoom+margin << endl
	  << "%%Magnification: 1.0000" << endl
	  << "%%EndComments" << endl
	  << "%%BeginProlog" << endl
	  << "%" << endl
	  << "/zoom " << zoom << " def" << endl
	  << "zoom zoom scale" << endl
	  << "/Helvetica findfont 7 zoom div scalefont setfont" << endl
	  << "0.06299 7.500 mul zoom div setlinewidth" << endl
	  << "/dashlen 4.5 zoom div def" << endl
	  << "/dashlen2 1.5 zoom div def" << endl
	  << "/ptlen 3.5 zoom div def" << endl << "%" << endl
	  << "% center-x center-y radius bp" << endl
	  << "% creates a path representing a circle" << endl
	  << "/bp { /radius exch def /cy exch def /cx exch def" << endl
	  << "radius 0 eq" << endl
	  << "{ newpath cx ptlen sub cy moveto ptlen ptlen add 0" << endl
	  << "rlineto cx cy ptlen sub moveto 0 ptlen ptlen add" << endl
	  << "rlineto }" << endl
	  << "{ newpath cx radius add cy moveto cx cy radius" << endl
	  << "0 360 arc closepath } ifelse } def" << "%" << endl
	  << "%center-x center-y radius ball" << endl
	  << "% draws a circle, unfilled, no dash" << endl
	  << "/ball { bp stroke } def" << endl
	  << "% center-x center-y radius bball" << endl
	  << "% draws a circle, filled, no dash" << endl
	  << "/bball { /radius exch def /cy exch def /cx exch def" << endl
	  << "cx cy radius radius 0 eq" << endl
	  << "{ gsave bp 0.6 setgray stroke" << endl
	  << "newpath cx cy moveto closepath 1 setlinecap" << endl
	  << "0.0 setgray stroke grestore }" << endl
	  << "{ bp gsave 0.6 setgray fill grestore stroke }" << endl
	  << "ifelse } def" << endl
	  << "% center-x center-y radius dball" << endl
	  << "% draws a circle, unfilled, dashed" << endl
	  << "/dball { gsave bp [dashlen] 0 setdash stroke grestore" << endl
	  << "} def" << endl
	  << "/eball { gsave bp 0.6 setgray stroke grestore } def" << endl
	  << "%" << endl
	  << "% center-x center-y radius label dist angle drawl" << endl
	  << "% Draws the label label." << endl
	  << "/drawl { /angle exch def /dist exch def /label exch def" << endl
	  << "/radius exch def /cy exch def /cx exch def" << endl
	  << "cx radius dist zoom div add angle cos mul add" << endl
	  << "cy radius dist zoom div add angle sin mul add" << endl
	  << "moveto label show } def" << endl
	  << "%" << endl
	  << "% Usage a b c find-roots x1 x2" << endl
	  << "% Finds the two roots of a x^2 + b x + c = 0, assuming" << endl
	  << "% that a is nonzero." << endl
	  << "/find-roots {" << endl
	  << "    % compute discriminant:" << endl
	  << "    3 copy 3 2 roll" << endl
	  << "    mul 4.0 mul neg exch dup mul add" << endl
	  << "    % stack: a b c disc" << endl
	  << "    sqrt" << endl