Lorenz.h 2.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109
  1. /*
  2. dsp/Lorenz.h
  3. Copyright 2001-4 Tim Goetze <tim@quitte.de>
  4. http://quitte.de/dsp/
  5. Lorenz fractal.
  6. */
  7. /*
  8. This program is free software; you can redistribute it and/or
  9. modify it under the terms of the GNU General Public License
  10. as published by the Free Software Foundation; either version 2
  11. of the License, or (at your option) any later version.
  12. This program is distributed in the hope that it will be useful,
  13. but WITHOUT ANY WARRANTY; without even the implied warranty of
  14. MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  15. GNU General Public License for more details.
  16. You should have received a copy of the GNU General Public License
  17. along with this program; if not, write to the Free Software
  18. Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
  19. 02111-1307, USA or point your web browser to http://www.gnu.org.
  20. */
  21. #ifndef _DSP_LORENZ_H_
  22. #define _DSP_LORENZ_H_
  23. namespace DSP {
  24. class Lorenz
  25. {
  26. public:
  27. double x[2], y[2], z[2];
  28. double h, a, b, c;
  29. int I;
  30. public:
  31. Lorenz()
  32. {
  33. h = 0.001;
  34. a = 10.0;
  35. b = 28.0;
  36. c = 8.0 / 3.0;
  37. }
  38. /* rate is normalized (0 .. 1) */
  39. void set_rate (double r)
  40. {
  41. h = max (.0000001, r * .015);
  42. }
  43. void init (double _h = .001, double seed = .0)
  44. {
  45. I = 0;
  46. x[0] = .1 + seed - frandom() * .1;
  47. y[0] = 0;
  48. z[0] = 0;
  49. /* progress quickly to get a 'stable' system */
  50. h = .001;
  51. int n = 10000 + min ((int) (10000 * seed), 10000);
  52. for (int i = 0; i < n; ++i)
  53. step();
  54. h = _h;
  55. }
  56. sample_t get()
  57. {
  58. step();
  59. return .5 * get_y() + get_z();
  60. }
  61. void step()
  62. {
  63. int J = I ^ 1;
  64. x[J] = x[I] + h * a * (y[I] - x[I]);
  65. y[J] = y[I] + h * (x[I] * (b - z[I]) - y[I]);
  66. z[J] = z[I] + h * (x[I] * y[I] - c * z[I]);
  67. I = J;
  68. }
  69. double get_x()
  70. {
  71. return .024 * (x[I] - .172);
  72. }
  73. double get_y()
  74. {
  75. return .018 * (y[I] - .172);
  76. }
  77. double get_z()
  78. {
  79. return .019 * (z[I] - 25.43);
  80. }
  81. };
  82. } /* namespace DSP */
  83. #endif /* _DSP_LORENZ_H_ */