src/tests/bezier.cpp

   1 /**
   2  * Tester
   3  * Bezier Curves - Compute them and put on Bitmaps
   4  */
   5 #include "main.h"
   6 #include <cmath>
   7 #include <ctime> // for performance measurements
   8 #include <unordered_map> // hashtable
   9 #include <map>
  10 #include "screen.h"
  11
  12 #include "pointstobitmap.h"
  13
  14 using namespace std;
  15 using namespace IPDF;
  16
  17 /**
  18  * Factorial
  19  * Use dynamic programming / recursion
  20  */
  21 int Factorial(int n)
  22 {
  23         static unordered_map<int, int> dp;
  24         static bool init = false;
  25         if (!init)
  26         {
  27                 init = true;
  28                 dp[0] = 1;
  29         }
  30         auto it = dp.find(n);
  31         if (it != dp.end())
  32                 return it->second;
  33         int result = n*Factorial(n-1);
  34         dp[n] = result;
  35         return result;
  36 }
  37
  38 /**
  39  * Binomial coefficients
  40  */
  41 int BinomialCoeff(int n, int k)
  42 {
  43         return Factorial(n) / Factorial(k) / Factorial(n-k);
  44 }
  45
  46 /**
  47  * Bernstein Basis Polynomial
  48  */
  49 template <class T> T Bernstein(int k, int n, const T & u)
  50 {
  51         return T(BinomialCoeff(n, k)) * pow(u, k) * pow(T(1.0) - u, n-k);
  52 }
  53
  54 /**
  55  * Evaluate a Bezier at a single point u
  56  */
  57 template <class T> pair<T, T> Bezier(const vector<pair<T, T> > & control, const T & u)
  58 {
  59         pair<T,T> result(0,0);
  60         for (size_t k = 0; k < control.size(); ++k)
  61         {
  62                 T bez(Bernstein<T>(k, control.size()-1, u));
  63                 result.first += control[k].first * bez;
  64                 result.second += control[k].second * bez;
  65         }
  66         return result;
  67 }
  68
  69 /**
  70  * Form a Bezier curve from the control points using the specified number of intervals
  71  */
  72 template <class T> vector<pair<T, T> > BezierCurve(const vector<pair<T, T> > & control, int intervals = 200)
  73 {
  74         vector<pair<T, T> > result(intervals+1);
  75         T dU = T(1)/intervals;
  76         T u = 0;
  77         for (int i = 0; i <= intervals; ++i)
  78         {
  79                 result[i] = Bezier<T>(control, u);
  80                 u += dU;
  81         }
  82         return result;
  83 }
  84
  85
  86
  87 /**
  88  * Make a circle using Beziers
  89  * Actually just makes one Bezier for a quadrant and then mirrors it for the others
  90  */
  91 template <class T> vector<pair<T, T> > BezierCircle(const T & scale = 1, int intervals = 50)
  92 {
  93         T k = T(4) * (pow(T(2), T(0.5)) - T(1)) / T(3); // k = 4/3 * (2^(1/2) - 1) because maths says so
  94         vector<pair<T, T> >control(4);
  95         control[0] = pair<T,T>(0,scale);
  96         control[1] = pair<T,T>(k*scale, scale);
  97         control[2] = pair<T,T>(scale, k*scale);
  98         control[3] = pair<T,T>(scale, 0);
  99
 100         auto points = BezierCurve<T>(control, intervals);
 101
 102         points.reserve(intervals*4);
 103         for (int i = 0; i < intervals; ++i)
 104         {
 105                 points.push_back(pair<T,T>(-points[i].first, points[i].second));
 106         }
 107         for (int i = 0; i < intervals; ++i)
 108         {
 109                 points.push_back(pair<T,T>(points[i].first, -points[i].second));
 110         }
 111         for (int i = 0; i < intervals; ++i)
 112         {
 113                 points.push_back(pair<T,T>(-points[i].first, -points[i].second));
 114         }
 115
 116         return points;
 117 }
 118
 119 /**
 120  * Test Bezier Curve with a scale applied to the *control* points
 121  */
 122 template <class T> vector<pair<T,T> > TestCurve(const T & scale = 1, long intervals = 50)
 123 {
 124         vector<pair<T,T> > control(3);
 125         control[0] = pair<T,T>(0,0);
 126         control[1] = pair<T,T>(0.4*scale, 0.8*scale);
 127         control[2] = pair<T,T>(scale,scale);
 128         return BezierCurve<T>(control, intervals);
 129 }
 130
 131
 132 /**
 133  * Test the Beziers over a range of scales. Can print the Nth bezier to stdout
 134  */
 135 void TestBeziers(long intervals = 50, long double start = 1e-38, long double end = 1e38, int print = -1)
 136 {
 137         int count = 0;
 138         // Choose our Bezier here
 139         //#define TESTBEZIER TestCurve
 140         #define TESTBEZIER BezierCircle
 141
 142         for (auto scale = start; scale < end; scale *= 10)
 143         {
 144                 auto f = TESTBEZIER<float>(scale, intervals);
 145                 auto d = TESTBEZIER<double>(scale, intervals);
 146                 auto l = TESTBEZIER<long double>(scale, intervals);
 147
 148                 // Make bitmap(s)
 149                 stringstream s;
 150                 s << "bezier" << count << "_" << scale << ".bmp";
 151                 PointsToBitmap<float>(f, scale, 5*(count++),0,0,s.str().c_str(),false);
 152
 153                 #if REAL > REAL_LONG_DOUBLE
 154                         auto r = TESTBEZIER<Real>(scale, intervals);
 155                 #endif
 156                 if (count++ == print)
 157                 {
 158                         for (auto i = 0; i < f.size(); ++i)
 159                         {
 160                                 printf("%d\t%.50lf\t%.50lf\t%.50llf\t%.50llf", f[i].first, f[i].second, d[i].first, d[i].second, l[i].first, l[i].second);
 161                                 #if REAL > REAL_LONG_DOUBLE
 162                                         printf("\t%.50lf\t%.50lf", r[i].first, r[i].second);
 163                                 #endif
 164                                 printf("\n");
 165                         }
 166                 }
 167
 168                 // optional: Overlay images instead
 169         }
 170 }
 171
 172 /**
 173  * main
 174  */
 175 int main(int argc, char ** argv)
 176 {
 177         TestBeziers(200);
 178         // Makes a movie, not really that exciting though
 179         //system("for i in *.bmp; do convert $i -filter Point -resize x600 +antialias $i; done");
 180         //system("ffmpeg -i bezier%d.bmp -framerate 2 -vcodec mpeg4 bezier.mp4");
 181 }
 182