synfig-core/src/modules/mod_noise/random_noise.cpp

   1 /* === S Y N F I G ========================================================= */
   2 /*!     \file mod_noise/random_noise.cpp
   3 **      \brief blehh
   4 **
   5 **      $Id$
   6 **
   7 **      \legal
   8 **      Copyright (c) 2002-2005 Robert B. Quattlebaum Jr., Adrian Bentley
   9 **      Copyright (c) 2007 Chris Moore
  10 **
  11 **      This package is free software; you can redistribute it and/or
  12 **      modify it under the terms of the GNU General Public License as
  13 **      published by the Free Software Foundation; either version 2 of
  14 **      the License, or (at your option) any later version.
  15 **
  16 **      This package is distributed in the hope that it will be useful,
  17 **      but WITHOUT ANY WARRANTY; without even the implied warranty of
  18 **      MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
  19 **      General Public License for more details.
  20 **      \endlegal
  21 */
  22 /* ========================================================================= */
  23
  24 /* === H E A D E R S ======================================================= */
  25
  26 #ifdef USING_PCH
  27 #       include "pch.h"
  28 #else
  29 #ifdef HAVE_CONFIG_H
  30 #       include <config.h>
  31 #endif
  32
  33 #include <synfig/general.h>
  34 #include "random_noise.h"
  35 #include <synfig/quick_rng.h>
  36 #include <cmath>
  37 #include <cstdlib>
  38 #endif
  39
  40 /* === M A C R O S ========================================================= */
  41
  42 #define PI      (3.1415927)
  43
  44 /* === G L O B A L S ======================================================= */
  45
  46 /* === P R O C E D U R E S ================================================= */
  47
  48 /* === M E T H O D S ======================================================= */
  49
  50 void
  51 RandomNoise::set_seed(int x)
  52 {
  53         seed_=x;
  54 }
  55
  56 float
  57 RandomNoise::operator()(const int salt,const int x,const int y,const int t)const
  58 {
  59         static const unsigned int a(21870);
  60         static const unsigned int b(11213);
  61         static const unsigned int c(36979);
  62         static const unsigned int d(31337);
  63
  64         quick_rng rng(
  65                 ( static_cast<unsigned int>(x+y)        * a ) ^
  66                 ( static_cast<unsigned int>(y+t)        * b ) ^
  67                 ( static_cast<unsigned int>(t+x)        * c ) ^
  68                 ( static_cast<unsigned int>(seed_+salt) * d )
  69         );
  70
  71         return rng.f() * 2.0f - 1.0f;
  72 }
  73
  74 float
  75 RandomNoise::operator()(SmoothType smooth,int subseed,float xf,float yf,float tf,int loop)const
  76 {
  77         int x((int)floor(xf));
  78         int y((int)floor(yf));
  79         int t((int)floor(tf));
  80         int t_1, t0, t1, t2;
  81
  82         if (loop)
  83         {
  84                 t0  = t % loop; if (t0  <  0   ) t0  += loop;
  85                 t_1 = t0 - 1;   if (t_1 <  0   ) t_1 += loop;
  86                 t1  = t0 + 1;   if (t1  >= loop) t1  -= loop;
  87                 t2  = t1 + 1;   if (t2  >= loop) t2  -= loop;
  88         }
  89         else
  90         {
  91                 t0  = t;
  92                 t_1 = t - 1;
  93                 t1  = t + 1;
  94                 t2  = t + 2;
  95         }
  96
  97         // synfig::info("%s:%d tf %.2f loop %d fraction %.2f ( -1,0,1,2 : %2d %2d %2d %2d)", __FILE__, __LINE__, tf, loop, tf-t, t_1, t0, t1, t2);
  98
  99         switch(smooth)
 100         {
 101         case SMOOTH_CUBIC:      // cubic
 102                 {
 103                         #define f(j,i,k)        ((*this)(subseed,i,j,k))
 104                         //Using catmull rom interpolation because it doesn't blur at all
 105                         // ( http://www.gamedev.net/reference/articles/article1497.asp )
 106                         //bezier curve with intermediate ctrl pts: 0.5/3(p(i+1) - p(i-1)) and similar
 107                         float xfa [4], tfa[4];
 108
 109                         //precalculate indices (all clamped) and offset
 110                         const int xa[] = {x-1,x,x+1,x+2};
 111                         const int ya[] = {y-1,y,y+1,y+2};
 112                         const int ta[] = {t_1,t0,t1,t2};
 113
 114                         const float dx(xf-x);
 115                         const float dy(yf-y);
 116                         const float dt(tf-t);
 117
 118                         //figure polynomials for each point
 119                         const float txf[] =
 120                         {
 121                                 0.5*dx*(dx*(dx*(-1) + 2) - 1),  //-t + 2t^2 -t^3
 122                                 0.5*(dx*(dx*(3*dx - 5)) + 2),   //2 - 5t^2 + 3t^3
 123                                 0.5*dx*(dx*(-3*dx + 4) + 1),    //t + 4t^2 - 3t^3
 124                                 0.5*dx*dx*(dx-1)                                //-t^2 + t^3
 125                         };
 126
 127                         const float tyf[] =
 128                         {
 129                                 0.5*dy*(dy*(dy*(-1) + 2) - 1),  //-t + 2t^2 -t^3
 130                                 0.5*(dy*(dy*(3*dy - 5)) + 2),   //2 - 5t^2 + 3t^3
 131                                 0.5*dy*(dy*(-3*dy + 4) + 1),    //t + 4t^2 - 3t^3
 132                                 0.5*dy*dy*(dy-1)                                //-t^2 + t^3
 133                         };
 134
 135                         const float ttf[] =
 136                         {
 137                                 0.5*dt*(dt*(dt*(-1) + 2) - 1),  //-t + 2t^2 -t^3
 138                                 0.5*(dt*(dt*(3*dt - 5)) + 2),   //2 - 5t^2 + 3t^3
 139                                 0.5*dt*(dt*(-3*dt + 4) + 1),    //t + 4t^2 - 3t^3
 140                                 0.5*dt*dt*(dt-1)                                //-t^2 + t^3
 141                         };
 142
 143                         //evaluate polynomial for each row
 144                         for(int i = 0; i < 4; ++i)
 145                         {
 146                                 for(int j = 0; j < 4; ++j)
 147                                 {
 148                                         tfa[j] = f(ya[i],xa[j],ta[0])*ttf[0] + f(ya[i],xa[j],ta[1])*ttf[1] + f(ya[i],xa[j],ta[2])*ttf[2] + f(ya[i],xa[j],ta[3])*ttf[3];
 149                                 }
 150                                 xfa[i] = tfa[0]*txf[0] + tfa[1]*txf[1] + tfa[2]*txf[2] + tfa[3]*txf[3];
 151                         }
 152
 153                         //return the cumulative column evaluation
 154                         return xfa[0]*tyf[0] + xfa[1]*tyf[1] + xfa[2]*tyf[2] + xfa[3]*tyf[3];
 155 #undef f
 156                 }
 157                 break;
 158
 159
 160         case SMOOTH_FAST_SPLINE:        // Fast Spline (non-animated)
 161                 {
 162 #define P(x)            (((x)>0)?((x)*(x)*(x)):0.0f)
 163 #define R(x)            ( P(x+2) - 4.0f*P(x+1) + 6.0f*P(x) - 4.0f*P(x-1) )*(1.0f/6.0f)
 164 #define F(i,j)          ((*this)(subseed,i+x,j+y)*(R((i)-a)*R(b-(j))))
 165 #define FT(i,j,k,l)     ((*this)(subseed,i+x,j+y,l)*(R((i)-a)*R(b-(j))*R((k)-c)))
 166 #define Z(i,j)          ret+=F(i,j)
 167 #define ZT(i,j,k,l) ret+=FT(i,j,k,l)
 168 #define X(i,j)          // placeholder... To make box more symmetric
 169 #define XT(i,j,k,l)     // placeholder... To make box more symmetric
 170
 171                 float a(xf-x), b(yf-y);
 172
 173                 // Interpolate
 174                 float ret(F(0,0));
 175                 Z(-1,-1); Z(-1, 0); Z(-1, 1); Z(-1, 2);
 176                 Z( 0,-1); X( 0, 0); Z( 0, 1); Z( 0, 2);
 177                 Z( 1,-1); Z( 1, 0); Z( 1, 1); Z( 1, 2);
 178                 Z( 2,-1); Z( 2, 0); Z( 2, 1); Z( 2, 2);
 179
 180                 return ret;
 181         }
 182
 183         case SMOOTH_SPLINE:     // Spline (animated)
 184                 {
 185                         float a(xf-x), b(yf-y), c(tf-t);
 186
 187                         // Interpolate
 188                         float ret(FT(0,0,0,t0));
 189                         ZT(-1,-1,-1,t_1); ZT(-1, 0,-1,t_1); ZT(-1, 1,-1,t_1); ZT(-1, 2,-1,t_1);
 190                         ZT( 0,-1,-1,t_1); ZT( 0, 0,-1,t_1); ZT( 0, 1,-1,t_1); ZT( 0, 2,-1,t_1);
 191                         ZT( 1,-1,-1,t_1); ZT( 1, 0,-1,t_1); ZT( 1, 1,-1,t_1); ZT( 1, 2,-1,t_1);
 192                         ZT( 2,-1,-1,t_1); ZT( 2, 0,-1,t_1); ZT( 2, 1,-1,t_1); ZT( 2, 2,-1,t_1);
 193
 194                         ZT(-1,-1, 0,t0 ); ZT(-1, 0, 0,t0 ); ZT(-1, 1, 0,t0 ); ZT(-1, 2, 0,t0 );
 195                         ZT( 0,-1, 0,t0 ); XT( 0, 0, 0,t0 ); ZT( 0, 1, 0,t0 ); ZT( 0, 2, 0,t0 );
 196                         ZT( 1,-1, 0,t0 ); ZT( 1, 0, 0,t0 ); ZT( 1, 1, 0,t0 ); ZT( 1, 2, 0,t0 );
 197                         ZT( 2,-1, 0,t0 ); ZT( 2, 0, 0,t0 ); ZT( 2, 1, 0,t0 ); ZT( 2, 2, 0,t0 );
 198
 199                         ZT(-1,-1, 1,t1 ); ZT(-1, 0, 1,t1 ); ZT(-1, 1, 1,t1 ); ZT(-1, 2, 1,t1 );
 200                         ZT( 0,-1, 1,t1 ); ZT( 0, 0, 1,t1 ); ZT( 0, 1, 1,t1 ); ZT( 0, 2, 1,t1 );
 201                         ZT( 1,-1, 1,t1 ); ZT( 1, 0, 1,t1 ); ZT( 1, 1, 1,t1 ); ZT( 1, 2, 1,t1 );
 202                         ZT( 2,-1, 1,t1 ); ZT( 2, 0, 1,t1 ); ZT( 2, 1, 1,t1 ); ZT( 2, 2, 1,t1 );
 203
 204                         ZT(-1,-1, 2,t2 ); ZT(-1, 0, 2,t2 ); ZT(-1, 1, 2,t2 ); ZT(-1, 2, 2,t2 );
 205                         ZT( 0,-1, 2,t2 ); ZT( 0, 0, 2,t2 ); ZT( 0, 1, 2,t2 ); ZT( 0, 2, 2,t2 );
 206                         ZT( 1,-1, 2,t2 ); ZT( 1, 0, 2,t2 ); ZT( 1, 1, 2,t2 ); ZT( 1, 2, 2,t2 );
 207                         ZT( 2,-1, 2,t2 ); ZT( 2, 0, 2,t2 ); ZT( 2, 1, 2,t2 ); ZT( 2, 2, 2,t2 );
 208
 209                         return ret;
 210
 211 /*
 212
 213                         float dx=xf-x;
 214                         float dy=yf-y;
 215                         float dt=tf-t;
 216
 217                         float ret=0;
 218                         int i,j,h;
 219                         for(h=-1;h<=2;h++)
 220                                 for(i=-1;i<=2;i++)
 221                                         for(j=-1;j<=2;j++)
 222                                                 ret+=(*this)(subseed,i+x,j+y,h+t)*(R(i-dx)*R(j-dy)*R(h-dt));
 223                         return ret;
 224 */
 225                 }
 226                 break;
 227 #undef X
 228 #undef Z
 229 #undef F
 230 #undef P
 231 #undef R
 232
 233         case SMOOTH_COSINE:
 234         if((float)t==tf)
 235         {
 236                 int x((int)floor(xf));
 237                 int y((int)floor(yf));
 238                 float a=xf-x;
 239                 float b=yf-y;
 240                 a=(1.0f-cos(a*PI))*0.5f;
 241                 b=(1.0f-cos(b*PI))*0.5f;
 242                 float c=1.0-a;
 243                 float d=1.0-b;
 244                 int x2=x+1,y2=y+1;
 245                 return
 246                         (*this)(subseed,x,y,t0)*(c*d)+
 247                         (*this)(subseed,x2,y,t0)*(a*d)+
 248                         (*this)(subseed,x,y2,t0)*(c*b)+
 249                         (*this)(subseed,x2,y2,t0)*(a*b);
 250         }
 251         else
 252         {
 253                 float a=xf-x;
 254                 float b=yf-y;
 255                 float c=tf-t;
 256
 257                 a=(1.0f-cos(a*PI))*0.5f;
 258                 b=(1.0f-cos(b*PI))*0.5f;
 259
 260                 // We don't perform this on the time axis, otherwise we won't
 261                 // get smooth motion
 262                 //c=(1.0f-cos(c*PI))*0.5f;
 263
 264                 float d=1.0-a;
 265                 float e=1.0-b;
 266                 float f=1.0-c;
 267
 268                 int x2=x+1,y2=y+1;
 269
 270                 return
 271                         (*this)(subseed,x,y,t0)*(d*e*f)+
 272                         (*this)(subseed,x2,y,t0)*(a*e*f)+
 273                         (*this)(subseed,x,y2,t0)*(d*b*f)+
 274                         (*this)(subseed,x2,y2,t0)*(a*b*f)+
 275                         (*this)(subseed,x,y,t1)*(d*e*c)+
 276                         (*this)(subseed,x2,y,t1)*(a*e*c)+
 277                         (*this)(subseed,x,y2,t1)*(d*b*c)+
 278                         (*this)(subseed,x2,y2,t1)*(a*b*c);
 279         }
 280         case SMOOTH_LINEAR:
 281         if((float)t==tf)
 282         {
 283                 int x((int)floor(xf));
 284                 int y((int)floor(yf));
 285                 float a=xf-x;
 286                 float b=yf-y;
 287                 float c=1.0-a;
 288                 float d=1.0-b;
 289                 int x2=x+1,y2=y+1;
 290                 return
 291                         (*this)(subseed,x,y,t0)*(c*d)+
 292                         (*this)(subseed,x2,y,t0)*(a*d)+
 293                         (*this)(subseed,x,y2,t0)*(c*b)+
 294                         (*this)(subseed,x2,y2,t0)*(a*b);
 295         }
 296         else
 297         {
 298
 299                 float a=xf-x;
 300                 float b=yf-y;
 301                 float c=tf-t;
 302
 303                 float d=1.0-a;
 304                 float e=1.0-b;
 305                 float f=1.0-c;
 306
 307                 int x2=x+1,y2=y+1;
 308
 309                 return
 310                         (*this)(subseed,x,y,t0)*(d*e*f)+
 311                         (*this)(subseed,x2,y,t0)*(a*e*f)+
 312                         (*this)(subseed,x,y2,t0)*(d*b*f)+
 313                         (*this)(subseed,x2,y2,t0)*(a*b*f)+
 314                         (*this)(subseed,x,y,t1)*(d*e*c)+
 315                         (*this)(subseed,x2,y,t1)*(a*e*c)+
 316                         (*this)(subseed,x,y2,t1)*(d*b*c)+
 317                         (*this)(subseed,x2,y2,t1)*(a*b*c);
 318         }
 319         default:
 320         case SMOOTH_DEFAULT:
 321                 return (*this)(subseed,x,y,t0);
 322         }
 323 }