1 /* === S Y N F I G ========================================================= */
5 ** $Id: random.cpp,v 1.6 2005/01/17 02:00:19 darco Exp $
8 ** Copyright (c) 2002 Robert B. Quattlebaum Jr.
10 ** This software and associated documentation
11 ** are CONFIDENTIAL and PROPRIETARY property of
12 ** the above-mentioned copyright holder.
14 ** You may not copy, print, publish, or in any
15 ** other way distribute this software without
16 ** a prior written agreement with
17 ** the copyright holder.
20 /* ========================================================================= */
22 /* === H E A D E R S ======================================================= */
36 // A fast 32-bit linear congruential random number generator
41 quick_rng(unsigned long seed):next(seed) { }
43 void set_seed(unsigned long x)
50 static const unsigned long a(1664525);
51 static const unsigned long c(1013904223);
63 static const float m(int(65535));
65 return float(i16())/m;
69 /* === M A C R O S ========================================================= */
71 #define PI (3.1415927)
73 /* === G L O B A L S ======================================================= */
75 /* === P R O C E D U R E S ================================================= */
77 /* === M E T H O D S ======================================================= */
80 Random::set_seed(int x)
86 Random::operator()(const int salt,const int x,const int y,const int t)const
88 static const unsigned int a(21870);
89 static const unsigned int b(11213);
90 static const unsigned int c(36979);
91 static const unsigned int d(31337);
94 ( static_cast<unsigned int>(x+y) * a ) ^
95 ( static_cast<unsigned int>(y+t) * b ) ^
96 ( static_cast<unsigned int>(t+x) * c ) ^
97 ( static_cast<unsigned int>(seed_+salt) * d )
100 return rng.f() * 2.0f - 1.0f;
104 Random::operator()(int smooth,int subseed,float xf,float yf,float tf)const
106 int x((int)floor(xf));
107 int y((int)floor(yf));
108 int t((int)floor(tf));
114 #define f(j,i,k) ((*this)(subseed,i,j,k))
115 //Using catmull rom interpolation because it doesn't blur at all
116 //bezier curve with intermediate ctrl pts: 0.5/3(p(i+1) - p(i-1)) and similar
117 float xfa [4], tfa[4];
119 //precalculate indices (all clamped) and offset
120 const int xa[] = {x-1,x,x+1,x+2};
122 const int ya[] = {y-1,y,y+1,y+2};
124 const int ta[] = {t-1,t,t+1,t+2};
126 const float dx(xf-x);
127 const float dy(yf-y);
128 const float dt(tf-t);
130 //figure polynomials for each point
133 0.5*dx*(dx*(dx*(-1) + 2) - 1), //-t + 2t^2 -t^3
134 0.5*(dx*(dx*(3*dx - 5)) + 2), //2 - 5t^2 + 3t^3
135 0.5*dx*(dx*(-3*dx + 4) + 1), //t + 4t^2 - 3t^3
136 0.5*dx*dx*(dx-1) //-t^2 + t^3
141 0.5*dy*(dy*(dy*(-1) + 2) - 1), //-t + 2t^2 -t^3
142 0.5*(dy*(dy*(3*dy - 5)) + 2), //2 - 5t^2 + 3t^3
143 0.5*dy*(dy*(-3*dy + 4) + 1), //t + 4t^2 - 3t^3
144 0.5*dy*dy*(dy-1) //-t^2 + t^3
149 0.5*dt*(dt*(dt*(-1) + 2) - 1), //-t + 2t^2 -t^3
150 0.5*(dt*(dt*(3*dt - 5)) + 2), //2 - 5t^2 + 3t^3
151 0.5*dt*(dt*(-3*dt + 4) + 1), //t + 4t^2 - 3t^3
152 0.5*dt*dt*(dt-1) //-t^2 + t^3
155 //evaluate polynomial for each row
156 for(int i = 0; i < 4; ++i)
158 for(int j = 0; j < 4; ++j)
160 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];
162 xfa[i] = tfa[0]*txf[0] + tfa[1]*txf[1] + tfa[2]*txf[2] + tfa[3]*txf[3];
165 //return the cumulative column evaluation
166 return xfa[0]*tyf[0] + xfa[1]*tyf[1] + xfa[2]*tyf[2] + xfa[3]*tyf[3];
172 case 5: // Fast Spline (non-animated)
174 #define P(x) (((x)>0)?((x)*(x)*(x)):0.0f)
175 #define R(x) ( P(x+2) - 4.0f*P(x+1) + 6.0f*P(x) - 4.0f*P(x-1) )*(1.0f/6.0f)
176 #define F(i,j) ((*this)(subseed,i+x,j+y)*(R((i)-a)*R(b-(j))))
177 #define FT(i,j,k) ((*this)(subseed,i+x,j+y,k+t)*(R((i)-a)*R(b-(j))*R((k)-c)))
178 #define Z(i,j) ret+=F(i,j)
179 #define ZT(i,j,k) ret+=FT(i,j,k)
180 #define X(i,j) // placeholder... To make box more symetric
181 #define XT(i,j,k) // placeholder... To make box more symetric
183 float a(xf-x), b(yf-y);
187 Z(-1,-1); Z(-1, 0); Z(-1, 1); Z(-1, 2);
188 Z( 0,-1); X( 0, 0); Z( 0, 1); Z( 0, 2);
189 Z( 1,-1); Z( 1, 0); Z( 1, 1); Z( 1, 2);
190 Z( 2,-1); Z( 2, 0); Z( 2, 1); Z( 2, 2);
195 case 3: // Spline (animated)
197 float a(xf-x), b(yf-y), c(tf-t);
200 float ret(FT(0,0,0));
201 ZT(-1,-1,-1); ZT(-1, 0,-1); ZT(-1, 1,-1); ZT(-1, 2,-1);
202 ZT( 0,-1,-1); ZT( 0, 0,-1); ZT( 0, 1,-1); ZT( 0, 2,-1);
203 ZT( 1,-1,-1); ZT( 1, 0,-1); ZT( 1, 1,-1); ZT( 1, 2,-1);
204 ZT( 2,-1,-1); ZT( 2, 0,-1); ZT( 2, 1,-1); ZT( 2, 2,-1);
206 ZT(-1,-1, 0); ZT(-1, 0, 0); ZT(-1, 1, 0); ZT(-1, 2, 0);
207 ZT( 0,-1, 0); XT( 0, 0, 0); ZT( 0, 1, 0); ZT( 0, 2, 0);
208 ZT( 1,-1, 0); ZT( 1, 0, 0); ZT( 1, 1, 0); ZT( 1, 2, 0);
209 ZT( 2,-1, 0); ZT( 2, 0, 0); ZT( 2, 1, 0); ZT( 2, 2, 0);
211 ZT(-1,-1, 1); ZT(-1, 0, 1); ZT(-1, 1, 1); ZT(-1, 2, 1);
212 ZT( 0,-1, 1); ZT( 0, 0, 1); ZT( 0, 1, 1); ZT( 0, 2, 1);
213 ZT( 1,-1, 1); ZT( 1, 0, 1); ZT( 1, 1, 1); ZT( 1, 2, 1);
214 ZT( 2,-1, 1); ZT( 2, 0, 1); ZT( 2, 1, 1); ZT( 2, 2, 1);
216 ZT(-1,-1, 2); ZT(-1, 0, 2); ZT(-1, 1, 2); ZT(-1, 2, 2);
217 ZT( 0,-1, 2); ZT( 0, 0, 2); ZT( 0, 1, 2); ZT( 0, 2, 2);
218 ZT( 1,-1, 2); ZT( 1, 0, 2); ZT( 1, 1, 2); ZT( 1, 2, 2);
219 ZT( 2,-1, 2); ZT( 2, 0, 2); ZT( 2, 1, 2); ZT( 2, 2, 2);
234 ret+=(*this)(subseed,i+x,j+y,h+t)*(R(i-dx)*R(j-dy)*R(h-dt));
248 int x((int)floor(xf));
249 int y((int)floor(yf));
252 a=(1.0f-cos(a*PI))*0.5f;
253 b=(1.0f-cos(b*PI))*0.5f;
258 (*this)(subseed,x,y,t)*(c*d)+
259 (*this)(subseed,x2,y,t)*(a*d)+
260 (*this)(subseed,x,y2,t)*(c*b)+
261 (*this)(subseed,x2,y2,t)*(a*b);
269 a=(1.0f-cos(a*3.1415927))*0.5f;
270 b=(1.0f-cos(b*3.1415927))*0.5f;
272 // We don't perform this on the time axis, otherwise we won't
274 //c=(1.0f-cos(c*3.1415927))*0.5f;
280 int x2=x+1,y2=y+1,t2=t+1;
283 (*this)(subseed,x,y,t)*(d*e*f)+
284 (*this)(subseed,x2,y,t)*(a*e*f)+
285 (*this)(subseed,x,y2,t)*(d*b*f)+
286 (*this)(subseed,x2,y2,t)*(a*b*f)+
287 (*this)(subseed,x,y,t2)*(d*e*c)+
288 (*this)(subseed,x2,y,t2)*(a*e*c)+
289 (*this)(subseed,x,y2,t2)*(d*b*c)+
290 (*this)(subseed,x2,y2,t2)*(a*b*c);
295 int x((int)floor(xf));
296 int y((int)floor(yf));
303 (*this)(subseed,x,y,t)*(c*d)+
304 (*this)(subseed,x2,y,t)*(a*d)+
305 (*this)(subseed,x,y2,t)*(c*b)+
306 (*this)(subseed,x2,y2,t)*(a*b);
319 int x2=x+1,y2=y+1,t2=t+1;
322 (*this)(subseed,x,y,t)*(d*e*f)+
323 (*this)(subseed,x2,y,t)*(a*e*f)+
324 (*this)(subseed,x,y2,t)*(d*b*f)+
325 (*this)(subseed,x2,y2,t)*(a*b*f)+
326 (*this)(subseed,x,y,t2)*(d*e*c)+
327 (*this)(subseed,x2,y,t2)*(a*e*c)+
328 (*this)(subseed,x,y2,t2)*(d*b*c)+
329 (*this)(subseed,x2,y2,t2)*(a*b*c);
333 return (*this)(subseed,x,y,t);