lemon/bits/bezier.h
author Peter Kovacs <kpeter@inf.elte.hu>
Fri, 03 Apr 2009 18:59:15 +0200
changeset 608 6ac5d9ae1d3d
parent 314 2cc60866a0c9
child 963 761fe0846f49
permissions -rw-r--r--
Support real types + numerical stability fix in NS (#254)

- Real types are supported by appropriate inicialization.
- A feature of the XTI spanning tree structure is removed to ensure
numerical stability (could cause problems using integer types).
The node potentials are updated always on the lower subtree,
in order to prevent overflow problems.
The former method isn't notably faster during to our tests.
alpar@209
     1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
alpar@128
     2
 *
alpar@209
     3
 * This file is a part of LEMON, a generic C++ optimization library.
alpar@128
     4
 *
alpar@440
     5
 * Copyright (C) 2003-2009
alpar@128
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
alpar@128
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
alpar@128
     8
 *
alpar@128
     9
 * Permission to use, modify and distribute this software is granted
alpar@128
    10
 * provided that this copyright notice appears in all copies. For
alpar@128
    11
 * precise terms see the accompanying LICENSE file.
alpar@128
    12
 *
alpar@128
    13
 * This software is provided "AS IS" with no warranty of any kind,
alpar@128
    14
 * express or implied, and with no claim as to its suitability for any
alpar@128
    15
 * purpose.
alpar@128
    16
 *
alpar@128
    17
 */
alpar@128
    18
alpar@128
    19
#ifndef LEMON_BEZIER_H
alpar@128
    20
#define LEMON_BEZIER_H
alpar@128
    21
kpeter@314
    22
//\ingroup misc
kpeter@314
    23
//\file
kpeter@314
    24
//\brief Classes to compute with Bezier curves.
kpeter@314
    25
//
kpeter@314
    26
//Up to now this file is used internally by \ref graph_to_eps.h
alpar@128
    27
alpar@128
    28
#include<lemon/dim2.h>
alpar@128
    29
alpar@128
    30
namespace lemon {
alpar@128
    31
  namespace dim2 {
alpar@128
    32
alpar@128
    33
class BezierBase {
alpar@128
    34
public:
alpar@184
    35
  typedef lemon::dim2::Point<double> Point;
alpar@128
    36
protected:
alpar@128
    37
  static Point conv(Point x,Point y,double t) {return (1-t)*x+t*y;}
alpar@128
    38
};
alpar@128
    39
alpar@128
    40
class Bezier1 : public BezierBase
alpar@128
    41
{
alpar@128
    42
public:
alpar@128
    43
  Point p1,p2;
alpar@128
    44
alpar@128
    45
  Bezier1() {}
alpar@128
    46
  Bezier1(Point _p1, Point _p2) :p1(_p1), p2(_p2) {}
alpar@209
    47
alpar@128
    48
  Point operator()(double t) const
alpar@128
    49
  {
alpar@128
    50
    //    return conv(conv(p1,p2,t),conv(p2,p3,t),t);
alpar@128
    51
    return conv(p1,p2,t);
alpar@128
    52
  }
alpar@128
    53
  Bezier1 before(double t) const
alpar@128
    54
  {
alpar@128
    55
    return Bezier1(p1,conv(p1,p2,t));
alpar@128
    56
  }
alpar@209
    57
alpar@128
    58
  Bezier1 after(double t) const
alpar@128
    59
  {
alpar@128
    60
    return Bezier1(conv(p1,p2,t),p2);
alpar@128
    61
  }
alpar@128
    62
alpar@128
    63
  Bezier1 revert() const { return Bezier1(p2,p1);}
alpar@128
    64
  Bezier1 operator()(double a,double b) const { return before(b).after(a/b); }
alpar@128
    65
  Point grad() const { return p2-p1; }
alpar@128
    66
  Point norm() const { return rot90(p2-p1); }
alpar@128
    67
  Point grad(double) const { return grad(); }
alpar@128
    68
  Point norm(double t) const { return rot90(grad(t)); }
alpar@128
    69
};
alpar@128
    70
alpar@128
    71
class Bezier2 : public BezierBase
alpar@128
    72
{
alpar@128
    73
public:
alpar@128
    74
  Point p1,p2,p3;
alpar@128
    75
alpar@128
    76
  Bezier2() {}
alpar@128
    77
  Bezier2(Point _p1, Point _p2, Point _p3) :p1(_p1), p2(_p2), p3(_p3) {}
alpar@128
    78
  Bezier2(const Bezier1 &b) : p1(b.p1), p2(conv(b.p1,b.p2,.5)), p3(b.p2) {}
alpar@128
    79
  Point operator()(double t) const
alpar@128
    80
  {
alpar@128
    81
    //    return conv(conv(p1,p2,t),conv(p2,p3,t),t);
alpar@128
    82
    return ((1-t)*(1-t))*p1+(2*(1-t)*t)*p2+(t*t)*p3;
alpar@128
    83
  }
alpar@128
    84
  Bezier2 before(double t) const
alpar@128
    85
  {
alpar@128
    86
    Point q(conv(p1,p2,t));
alpar@128
    87
    Point r(conv(p2,p3,t));
alpar@128
    88
    return Bezier2(p1,q,conv(q,r,t));
alpar@128
    89
  }
alpar@209
    90
alpar@128
    91
  Bezier2 after(double t) const
alpar@128
    92
  {
alpar@128
    93
    Point q(conv(p1,p2,t));
alpar@128
    94
    Point r(conv(p2,p3,t));
alpar@128
    95
    return Bezier2(conv(q,r,t),r,p3);
alpar@128
    96
  }
alpar@128
    97
  Bezier2 revert() const { return Bezier2(p3,p2,p1);}
alpar@128
    98
  Bezier2 operator()(double a,double b) const { return before(b).after(a/b); }
alpar@128
    99
  Bezier1 grad() const { return Bezier1(2.0*(p2-p1),2.0*(p3-p2)); }
alpar@128
   100
  Bezier1 norm() const { return Bezier1(2.0*rot90(p2-p1),2.0*rot90(p3-p2)); }
alpar@128
   101
  Point grad(double t) const { return grad()(t); }
alpar@128
   102
  Point norm(double t) const { return rot90(grad(t)); }
alpar@128
   103
};
alpar@128
   104
alpar@128
   105
class Bezier3 : public BezierBase
alpar@128
   106
{
alpar@128
   107
public:
alpar@128
   108
  Point p1,p2,p3,p4;
alpar@128
   109
alpar@128
   110
  Bezier3() {}
alpar@128
   111
  Bezier3(Point _p1, Point _p2, Point _p3, Point _p4)
alpar@128
   112
    : p1(_p1), p2(_p2), p3(_p3), p4(_p4) {}
alpar@209
   113
  Bezier3(const Bezier1 &b) : p1(b.p1), p2(conv(b.p1,b.p2,1.0/3.0)),
alpar@209
   114
                              p3(conv(b.p1,b.p2,2.0/3.0)), p4(b.p2) {}
alpar@128
   115
  Bezier3(const Bezier2 &b) : p1(b.p1), p2(conv(b.p1,b.p2,2.0/3.0)),
alpar@209
   116
                              p3(conv(b.p2,b.p3,1.0/3.0)), p4(b.p3) {}
alpar@209
   117
alpar@209
   118
  Point operator()(double t) const
alpar@128
   119
    {
alpar@128
   120
      //    return Bezier2(conv(p1,p2,t),conv(p2,p3,t),conv(p3,p4,t))(t);
alpar@128
   121
      return ((1-t)*(1-t)*(1-t))*p1+(3*t*(1-t)*(1-t))*p2+
alpar@209
   122
        (3*t*t*(1-t))*p3+(t*t*t)*p4;
alpar@128
   123
    }
alpar@128
   124
  Bezier3 before(double t) const
alpar@128
   125
    {
alpar@128
   126
      Point p(conv(p1,p2,t));
alpar@128
   127
      Point q(conv(p2,p3,t));
alpar@128
   128
      Point r(conv(p3,p4,t));
alpar@128
   129
      Point a(conv(p,q,t));
alpar@128
   130
      Point b(conv(q,r,t));
alpar@128
   131
      Point c(conv(a,b,t));
alpar@128
   132
      return Bezier3(p1,p,a,c);
alpar@128
   133
    }
alpar@209
   134
alpar@128
   135
  Bezier3 after(double t) const
alpar@128
   136
    {
alpar@128
   137
      Point p(conv(p1,p2,t));
alpar@128
   138
      Point q(conv(p2,p3,t));
alpar@128
   139
      Point r(conv(p3,p4,t));
alpar@128
   140
      Point a(conv(p,q,t));
alpar@128
   141
      Point b(conv(q,r,t));
alpar@128
   142
      Point c(conv(a,b,t));
alpar@128
   143
      return Bezier3(c,b,r,p4);
alpar@128
   144
    }
alpar@128
   145
  Bezier3 revert() const { return Bezier3(p4,p3,p2,p1);}
alpar@128
   146
  Bezier3 operator()(double a,double b) const { return before(b).after(a/b); }
alpar@128
   147
  Bezier2 grad() const { return Bezier2(3.0*(p2-p1),3.0*(p3-p2),3.0*(p4-p3)); }
alpar@128
   148
  Bezier2 norm() const { return Bezier2(3.0*rot90(p2-p1),
alpar@209
   149
                                  3.0*rot90(p3-p2),
alpar@209
   150
                                  3.0*rot90(p4-p3)); }
alpar@128
   151
  Point grad(double t) const { return grad()(t); }
alpar@128
   152
  Point norm(double t) const { return rot90(grad(t)); }
alpar@128
   153
alpar@128
   154
  template<class R,class F,class S,class D>
alpar@209
   155
  R recSplit(F &_f,const S &_s,D _d) const
alpar@128
   156
  {
alpar@128
   157
    const Point a=(p1+p2)/2;
alpar@128
   158
    const Point b=(p2+p3)/2;
alpar@128
   159
    const Point c=(p3+p4)/2;
alpar@128
   160
    const Point d=(a+b)/2;
alpar@128
   161
    const Point e=(b+c)/2;
alpar@128
   162
    const Point f=(d+e)/2;
alpar@128
   163
    R f1=_f(Bezier3(p1,a,d,e),_d);
alpar@128
   164
    R f2=_f(Bezier3(e,d,c,p4),_d);
alpar@128
   165
    return _s(f1,f2);
alpar@128
   166
  }
alpar@209
   167
alpar@128
   168
};
alpar@128
   169
alpar@128
   170
alpar@128
   171
} //END OF NAMESPACE dim2
alpar@128
   172
} //END OF NAMESPACE lemon
alpar@128
   173
alpar@128
   174
#endif // LEMON_BEZIER_H