1 | /* |
---|
2 | ----------------------------------------------------------------------------- |
---|
3 | This source file is part of OGRE |
---|
4 | (Object-oriented Graphics Rendering Engine) |
---|
5 | For the latest info, see http://www.ogre3d.org/ |
---|
6 | |
---|
7 | Copyright (c) 2000-2006 Torus Knot Software Ltd |
---|
8 | Also see acknowledgements in Readme.html |
---|
9 | |
---|
10 | You may use this sample code for anything you like, it is not covered by the |
---|
11 | LGPL like the rest of the engine. |
---|
12 | ----------------------------------------------------------------------------- |
---|
13 | */ |
---|
14 | /* Original author: John C. Hare |
---|
15 | * Date: June 26, 1996 (although I've had these around for at least 6 years) |
---|
16 | * Adapted to C++ by W.J. van der Laan 2004 |
---|
17 | */ |
---|
18 | #ifndef H_Q_Julia |
---|
19 | #define H_Q_Julia |
---|
20 | #include <cmath> |
---|
21 | |
---|
22 | /** |
---|
23 | * Simple, fast, inline quaternion math functions |
---|
24 | */ |
---|
25 | struct Quat { |
---|
26 | float r,i,j,k; |
---|
27 | }; |
---|
28 | |
---|
29 | inline void qadd(Quat &a, const Quat &b) |
---|
30 | { |
---|
31 | a.r += b.r; |
---|
32 | a.i += b.i; |
---|
33 | a.j += b.j; |
---|
34 | a.k += b.k; |
---|
35 | } |
---|
36 | |
---|
37 | inline void qmult(Quat &c, const Quat &a, const Quat &b) |
---|
38 | { |
---|
39 | c.r = a.r*b.r - a.i*b.i - a.j*b.j - a.k*b.k; |
---|
40 | c.i = a.r*b.i + a.i*b.r + a.j*b.k - a.k*b.j; |
---|
41 | c.j = a.r*b.j + a.j*b.r + a.k*b.i - a.i*b.k; |
---|
42 | c.k = a.r*b.k + a.k*b.r + a.i*b.j - a.j*b.i; |
---|
43 | } |
---|
44 | |
---|
45 | inline void qsqr(Quat &b, const Quat &a) |
---|
46 | { |
---|
47 | b.r = a.r*a.r - a.i*a.i - a.j*a.j - a.k*a.k; |
---|
48 | b.i = 2.0f*a.r*a.i; |
---|
49 | b.j = 2.0f*a.r*a.j; |
---|
50 | b.k = 2.0f*a.r*a.k; |
---|
51 | } |
---|
52 | |
---|
53 | /** |
---|
54 | * Implicit function that evaluates the Julia set. |
---|
55 | */ |
---|
56 | class Julia { |
---|
57 | private: |
---|
58 | float global_real, global_imag, global_theta; |
---|
59 | Quat oc,c,eio,emio; |
---|
60 | public: |
---|
61 | Julia(float global_real, float global_imag, float global_theta); |
---|
62 | inline float eval(float x, float y, float z) { |
---|
63 | Quat q, temp; |
---|
64 | int i; |
---|
65 | |
---|
66 | q.r = x; |
---|
67 | q.i = y; |
---|
68 | q.j = z; |
---|
69 | q.k = 0.0; |
---|
70 | |
---|
71 | for (i = 30; i > 0; i--) { |
---|
72 | qsqr(temp, q); |
---|
73 | qmult(q, emio, temp); |
---|
74 | qadd(q, c); |
---|
75 | |
---|
76 | if (q.r*q.r + q.i*q.i + q.j*q.j + q.k*q.k > 8.0) |
---|
77 | break; |
---|
78 | } |
---|
79 | |
---|
80 | return((float)i); |
---|
81 | } |
---|
82 | }; |
---|
83 | |
---|
84 | |
---|
85 | Julia::Julia(float global_real, float global_imag, float global_theta): |
---|
86 | global_real(global_real), global_imag(global_imag), global_theta(global_theta) { |
---|
87 | |
---|
88 | oc.r = global_real; |
---|
89 | oc.i = global_imag; |
---|
90 | oc.j = oc.k = 0.0; |
---|
91 | |
---|
92 | eio.r = cos(global_theta); |
---|
93 | eio.i = sin(global_theta); |
---|
94 | eio.j = 0.0; |
---|
95 | eio.k = 0.0; |
---|
96 | |
---|
97 | emio.r = cos(-global_theta); |
---|
98 | emio.i = sin(-global_theta); |
---|
99 | emio.j = 0.0; |
---|
100 | emio.k = 0.0; |
---|
101 | |
---|
102 | /*** |
---|
103 | *** multiply eio*c only once at the beginning of iteration |
---|
104 | *** (since q |-> sqrt(eio*(q-eio*c))) |
---|
105 | *** q -> e-io*q^2 - eio*c |
---|
106 | ***/ |
---|
107 | |
---|
108 | qmult(c, eio,oc); |
---|
109 | |
---|
110 | } |
---|
111 | |
---|
112 | |
---|
113 | #endif |
---|