git.mutantstargoat.com Git - summerhack/blob - tools/curve_draw/curves.cc

   1 /*
   2 Copyright 2004 John Tsiombikas <nuclear@siggraph.org>
   3
   4 This file is part of the graphics core library.
   5
   6 the graphics core library is free software; you can redistribute it and/or modify
   7 it under the terms of the GNU General Public License as published by
   8 the Free Software Foundation; either version 2 of the License, or
   9 (at your option) any later version.
  10
  11 the graphics core library is distributed in the hope that it will be useful,
  12 but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14 GNU General Public License for more details.
  15
  16 You should have received a copy of the GNU General Public License
  17 along with the graphics core library; if not, write to the Free Software
  18 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  19 */
  20
  21 /* 3D Curves
  22  *
  23  * author: John Tsiombikas 2003
  24  * modified:
  25  *              John Tsiombikas 2004
  26  *              Mihalis Georgoulopoulos 2004
  27  */
  28
  29 #include <cstdio>
  30 #include <cstring>
  31 #include <cmath>
  32 #include <cctype>
  33 #include <cassert>
  34 #include "curves.h"
  35 //#include "common/err_msg.h"
  36
  37 Curve::Curve() {
  38         arc_parametrize = false;
  39         ease_curve = 0;
  40         samples = 0;
  41
  42         set_ease_sample_count(100);
  43 }
  44
  45 Curve::~Curve() {
  46         delete [] samples;
  47
  48 }
  49
  50 void Curve::set_arc_parametrization(bool state) {
  51         arc_parametrize = state;
  52 }
  53
  54 #define Param   0
  55 #define ArcLen  1
  56
  57 void Curve::sample_arc_lengths() {
  58         const int samplesPerSegment = 30;
  59         sample_count = get_segment_count() * samplesPerSegment;
  60
  61         arc_parametrize = false;        // to be able to interpolate with the original values
  62
  63         samples = new Vector2[sample_count];
  64         Vector3 prevpos;
  65         scalar_t step = 1.0f / (scalar_t)(sample_count-1);
  66         for(int i=0; i<sample_count; i++) {
  67                 scalar_t t = step * (scalar_t)i;
  68                 Vector3 pos = interpolate(t);
  69                 samples[i][Param] = t;
  70                 if(!i) {
  71                         samples[i][ArcLen] = 0.0f;
  72                 } else {
  73                         samples[i][ArcLen] = (pos - prevpos).length() + samples[i-1][ArcLen];
  74                 }
  75                 prevpos = pos;
  76         }
  77
  78         // normalize arc lenghts
  79         scalar_t maxlen = samples[sample_count-1][ArcLen];
  80         for(int i=0; i<sample_count; i++) {
  81                 samples[i][ArcLen] /= maxlen;
  82         }
  83
  84         arc_parametrize = true;
  85 }
  86
  87 static int binary_search(Vector2 *array, scalar_t key, int begin, int end) {
  88         int middle = begin + ((end - begin)>>1);
  89
  90         if(array[middle][ArcLen] == key) return middle;
  91         if(end == begin) return middle;
  92
  93         if(key < array[middle][ArcLen]) return binary_search(array, key, begin, middle);
  94         if(key > array[middle][ArcLen]) return binary_search(array, key, middle+1, end);
  95         return -1;      // just to make the compiler shut the fuck up
  96 }
  97
  98 scalar_t Curve::parametrize(scalar_t t) const {
  99         if(!samples) const_cast<Curve*>(this)->sample_arc_lengths();
 100
 101         int samplepos = binary_search(samples, t, 0, sample_count);
 102         scalar_t par = samples[samplepos][Param];
 103         scalar_t len = samples[samplepos][ArcLen];
 104
 105         // XXX: I can't remember the significance of this condition, I had xsmall_number here
 106         // previously and if t was 0.9999 it broke. I just changed the number blindly which fixed
 107         // the breakage but I should investigate further at some point.
 108         if((len - t) < 0.0005) return par;
 109
 110         if(len < t) {
 111                 if(!samplepos) return par;
 112                 scalar_t prevlen = samples[samplepos-1][ArcLen];
 113                 scalar_t prevpar = samples[samplepos-1][Param];
 114                 scalar_t p = (t - prevlen) / (len - prevlen);
 115                 return prevpar + (par - prevpar) * p;
 116         } else {
 117                 if(samplepos >= sample_count) return par;
 118                 scalar_t nextlen = samples[samplepos+1][ArcLen];
 119                 scalar_t nextpar = samples[samplepos+1][Param];
 120                 scalar_t p = (t - len) / (nextlen - len);
 121                 return par + (nextpar - par) * p;
 122         }
 123
 124         return par;     // not gonna happen
 125 }
 126
 127
 128 #define MIN(a, b) ((a) < (b) ? (a) : (b))
 129 #define MAX(a, b) ((a) > (b) ? (a) : (b))
 130
 131 scalar_t Curve::ease(scalar_t t) const {
 132         if(!ease_curve) return t;
 133
 134         const_cast<Curve*>(this)->ease_curve->set_arc_parametrization(true);
 135         scalar_t et = ease_curve->interpolate(t).y;
 136
 137         return MIN(MAX(et, 0.0f), 1.0f);
 138 }
 139
 140
 141 void Curve::add_control_point(const Vector3 &cp) {
 142         control_points.push_back(cp);
 143         delete [] samples;
 144         samples = 0;
 145 }
 146
 147 void Curve::remove_control_point(int index) {
 148         if(index < 0 || index >= control_points.size()) return;
 149
 150         ListNode<Vector3> *node = control_points.begin();
 151         for(int i=0; i<index; i++) {
 152                 node = node->next;
 153         }
 154
 155         control_points.erase(node);
 156 }
 157
 158 #define MIN(a, b)       ((a) < (b) ? (a) : (b))
 159 #define MAX(a, b)       ((a) > (b) ? (a) : (b))
 160
 161 Vector3 *Curve::get_control_point(int index) {
 162         index = MAX(0, MIN(index, control_points.size() - 1));
 163
 164         ListNode<Vector3> *node = control_points.begin();
 165         for(int i=0; i<index; i++) {
 166                 assert(node);
 167                 node = node->next;
 168         }
 169
 170         assert(node);
 171         return &node->data;
 172 }
 173
 174 int Curve::get_point_count() const {
 175         return control_points.size();
 176 }
 177
 178 void Curve::set_ease_curve(Curve *curve) {
 179         ease_curve = curve;
 180 }
 181
 182 void Curve::set_ease_sample_count(int count) {
 183         ease_sample_count = count;
 184         ease_step = (int)(1.0f / (scalar_t)ease_sample_count);
 185 }
 186
 187 Vector3 Curve::operator ()(scalar_t t) const {
 188         return interpolate(t);
 189 }
 190
 191 ///////////////// B-Spline implementation ////////////////////
 192
 193 int BSplineCurve::get_segment_count() const {
 194         return control_points.size() - 3;
 195 }
 196
 197 Vector3 BSplineCurve::interpolate(scalar_t t) const {
 198         if(t > 1.0) t = 1.0;
 199         if(t < 0.0) t = 0.0;
 200
 201         if(control_points.size() < 4) return Vector3(0, 0, 0);
 202
 203         if(arc_parametrize) {
 204                 t = ease(parametrize(t));
 205         }
 206
 207         // find the appropriate segment of the spline that t lies and calculate the piecewise parameter
 208         t = (scalar_t)(control_points.size() - 3) * t;
 209         int seg = (int)t;
 210         t -= (scalar_t)floor(t);
 211         if(seg >= get_segment_count()) {
 212                 seg = get_segment_count() - 1;
 213                 t = 1.0f;
 214         }
 215
 216         ListNode<Vector3> *iter = const_cast<BSplineCurve*>(this)->control_points.begin();
 217         for(int i=0; i<seg; i++) iter = iter->next;
 218
 219         Vector3 Cp[4];
 220         for(int i=0; i<4; i++) {
 221         Cp[i] = iter->data;
 222                 iter = iter->next;
 223         }
 224
 225         Vector3 res;
 226         res.x = bspline(Cp[0].x, Cp[1].x, Cp[2].x, Cp[3].x, t);
 227         res.y = bspline(Cp[0].y, Cp[1].y, Cp[2].y, Cp[3].y, t);
 228         res.z = bspline(Cp[0].z, Cp[1].z, Cp[2].z, Cp[3].z, t);
 229
 230         return res;
 231 }
 232
 233 //////////////// Catmull-Rom Spline implementation //////////////////
 234
 235 int CatmullRomSplineCurve::get_segment_count() const {
 236         return control_points.size() - 1;
 237 }
 238
 239 Vector3 CatmullRomSplineCurve::interpolate(scalar_t t) const {
 240         if(t > 1.0) t = 1.0;
 241         if(t < 0.0) t = 0.0;
 242
 243         if(control_points.size() < 2) return Vector3(0, 0, 0);
 244
 245         if(arc_parametrize) {
 246                 t = ease(parametrize(t));
 247         }
 248
 249         // find the appropriate segment of the spline that t lies and calculate the piecewise parameter
 250         t = (scalar_t)(control_points.size() - 1) * t;
 251         int seg = (int)t;
 252         t -= (scalar_t)floor(t);
 253         if(seg >= get_segment_count()) {
 254                 seg = get_segment_count() - 1;
 255                 t = 1.0f;
 256         }
 257
 258         Vector3 Cp[4];
 259         ListNode<Vector3> *iter = const_cast<CatmullRomSplineCurve*>(this)->control_points.begin();
 260         for(int i=0; i<seg; i++) iter = iter->next;
 261
 262         Cp[1] = iter->data;
 263         Cp[2] = iter->next->data;
 264
 265         if(!seg) {
 266                 Cp[0] = Cp[1];
 267         } else {
 268                 Cp[0] = iter->prev->data;
 269         }
 270
 271         if(seg == control_points.size() - 2) {
 272                 Cp[3] = Cp[2];
 273         } else {
 274                 Cp[3] = iter->next->next->data;
 275         }
 276
 277         Vector3 res;
 278         res.x = catmull_rom_spline(Cp[0].x, Cp[1].x, Cp[2].x, Cp[3].x, t);
 279         res.y = catmull_rom_spline(Cp[0].y, Cp[1].y, Cp[2].y, Cp[3].y, t);
 280         res.z = catmull_rom_spline(Cp[0].z, Cp[1].z, Cp[2].z, Cp[3].z, t);
 281
 282         return res;
 283 }
 284
 285 /* BezierSpline implementation - (MG)
 286  */
 287 int BezierSpline::get_segment_count() const
 288 {
 289         if (control_points.size() < 0) return 0;
 290         return control_points.size() / 4;
 291 }
 292
 293 Vector3 BezierSpline::interpolate(scalar_t t) const
 294 {
 295         if (!get_segment_count()) return Vector3(0, 0, 0);
 296
 297         if (arc_parametrize)
 298         {
 299                 t = ease(parametrize(t));
 300         }
 301
 302         t = (t * get_segment_count());
 303         int seg = (int) t;
 304         t -= (scalar_t) floor(t);
 305         if (seg >= get_segment_count())
 306         {
 307                 seg = get_segment_count() - 1;
 308                 t = 1.0f;
 309         }
 310
 311         seg *= 4;
 312         ListNode<Vector3> *iter = const_cast<BezierSpline*>(this)->control_points.begin();
 313         for (int i = 0; i < seg; i++) iter = iter->next;
 314
 315         Vector3 Cp[4];
 316         for (int i = 0; i < 4; i++)
 317         {
 318                 Cp[i] = iter->data;
 319                 iter = iter->next;
 320         }
 321
 322         // interpolate
 323         return bezier(Cp[seg], Cp[seg + 1], Cp[seg + 2], Cp[seg + 3], t);
 324 }
 325
 326 Vector3 BezierSpline::get_tangent(scalar_t t)
 327 {
 328         if (!get_segment_count()) return Vector3(0, 0, 0);
 329
 330         if (arc_parametrize)
 331         {
 332                 t = ease(parametrize(t));
 333         }
 334
 335         t = (t * get_segment_count());
 336         int seg = (int) t;
 337         t -= (scalar_t) floor(t);
 338         if (seg >= get_segment_count())
 339         {
 340                 seg = get_segment_count() - 1;
 341                 t = 1.0f;
 342         }
 343
 344         seg *= 4;
 345         ListNode<Vector3> *iter = const_cast<BezierSpline*>(this)->control_points.begin();
 346         for (int i = 0; i < seg; i++) iter = iter->next;
 347
 348         Vector3 Cp[4];
 349         for (int i = 0; i < 4; i++)
 350         {
 351                 Cp[i] = iter->data;
 352                 iter = iter->next;
 353         }
 354
 355         // interpolate
 356         return bezier_tangent(Cp[0], Cp[1], Cp[2], Cp[3], t);
 357 }
 358
 359 Vector3 BezierSpline::get_control_point(int i) const
 360 {
 361
 362         ListNode<Vector3> *iter = const_cast<BezierSpline*>(this)->control_points.begin();
 363         for (int j = 0; j < i; j++)
 364         {
 365                 if (!iter->next)
 366                 {
 367                         return Vector3(0, 0, 0);
 368                 }
 369                 iter = iter->next;
 370         }
 371
 372         return iter->data;
 373 }
 374
 375
 376
 377 bool save_curve(const char *fname, const Curve *curve) {
 378         FILE *fp = fopen(fname, "w");
 379         if(!fp) {
 380                 fprintf(stderr, "failed to save the curve %s", curve->name.c_str());
 381                 return false;
 382         }
 383
 384         fputs("curve_3dengfx\n", fp);
 385         fputs(curve->name.c_str(), fp);
 386         fputs("\n", fp);
 387
 388         if(dynamic_cast<const BSplineCurve*>(curve)) {
 389                 fputs("bspline\n", fp);
 390         } else if(dynamic_cast<const CatmullRomSplineCurve*>(curve)) {
 391                 fputs("catmullrom\n", fp);
 392         } else if(dynamic_cast<const BezierSpline*>(curve)) {
 393                 fputs("bezier\n", fp);
 394         } else {
 395                 fprintf(stderr, "unknown spline type, save failed %s", curve->name.c_str());
 396                 fclose(fp);
 397                 remove(fname);
 398                 return false;
 399         }
 400
 401         fprintf(fp, "%d\n", curve->control_points.size());
 402
 403         const ListNode<Vector3> *node = curve->control_points.begin();
 404         while(node) {
 405                 fprintf(fp, "%f %f %f\n", node->data.x, node->data.y, node->data.z);
 406                 node = node->next;
 407         }
 408
 409         fclose(fp);
 410         return true;
 411 }
 412
 413 Curve *load_curve(const char *fname) {
 414         FILE *fp = fopen(fname, "r");
 415         if(!fp) {
 416                 fprintf(stderr, "failed to open file %s", fname);
 417                 return 0;
 418         }
 419
 420         char buffer[256];
 421
 422         fgets(buffer, 256, fp);
 423         if(strcmp(buffer, "curve_3dengfx\n") != 0) {
 424                 fprintf(stderr, "load_curve failed, %s is not a curve file", fname);
 425                 fclose(fp);
 426                 return 0;
 427         }
 428
 429         Curve *curve;
 430
 431         fgets(buffer, 256, fp);
 432         std::string name = buffer;
 433
 434         fgets(buffer, 256, fp);
 435         if(!strcmp(buffer, "bspline\n")) {
 436                 curve = new BSplineCurve;
 437         } else if(!strcmp(buffer, "catmullrom\n")) {
 438                 curve = new CatmullRomSplineCurve;
 439         } else /*if(!strcmp(buffer, "bezier"))*/ {
 440                 fprintf(stderr, "unsupported curve type (%s) or not a curve file", buffer);
 441                 fclose(fp);
 442                 return 0;
 443         }
 444
 445         curve->name = name;
 446
 447         fgets(buffer, 256, fp);
 448         if(!isdigit(buffer[0])) {
 449                 fprintf(stderr, "load_curve failed, %s is not a valid curve file (count: %s)", fname, buffer);
 450                 delete curve;
 451                 fclose(fp);
 452                 return 0;
 453         }
 454         int count = atoi(buffer);
 455
 456         int failed = count;
 457         for(int i=0; i<count; i++, failed--) {
 458                 fgets(buffer, 256, fp);
 459                 if(!isdigit(buffer[0]) && buffer[0] != '.' && buffer[0] != '-') {
 460                         break;
 461                 }
 462                 float x = atof(buffer);
 463
 464                 char *ptr = strchr(buffer, ' ');
 465                 if(!ptr || (!isdigit(ptr[1]) && ptr[1] != '.' && ptr[1] != '-')) {
 466                         break;
 467                 }
 468                 float y = atof(++ptr);
 469
 470                 ptr = strchr(ptr, ' ');
 471                 if(!ptr || (!isdigit(ptr[1]) && ptr[1] != '.' && ptr[1] != '-')) {
 472                         break;
 473                 }
 474                 float z = atof(++ptr);
 475
 476                 curve->add_control_point(Vector3(x, y, z));
 477         }
 478
 479         fclose(fp);
 480
 481         if(failed) {
 482                 fprintf(stderr, "load_curve failed to read the data, %s is not a valid curve file", fname);
 483                 delete curve;
 484                 return 0;
 485         }
 486
 487         return curve;
 488 }