}

Delphi.int.ru — Портал программистов

Вход Регистрация | Забыли пароль?

Просмотр кода

Идентификатор: 2b7d3f01 Описание: FFT Delphi Код загружен: 23 июня 2011, 18:03 (mirt.steelwater)

  1. type
  2. { массив действительных чисел }
  3.   PDoubleArray =
  4. PDoubleArray = ^TDoubleArray;
  5. TDoubleArray = array [0..0] of Double;
  6. {nction = ^TDoubleFunction;
  7.   TDoubleFunction = function (const anArray: PDoubleArray;
  8.   const X: WORD) : Double;
  9.   { параметры вызова функции }
  10.   PDoubleFunctionParams = ^TDoubleFun
  11. PDoubleFunctionParams = ^TDoubleFunctionParams;
  12. TDoubleFunctionParams = packed record
  13. X1 : WORD;
  14. X2 : WORD;
  15. dX : WORD;
  16. end;
  17.  
  18. function spectr_pool (const anArray: PDoubleArray;
  19. const X: WORD) : Double;
  20. begin
  21. Result := anArray^ [2*x];
  22. end;
  23.  
  24. function fft_pool (const anArray: PDoubleArray;
  25. const X: WORD) : Double;
  26. begin
  27. Result := sqr ( anArray^ [x] );
  28. end;
  29.  
  30. function entropy_pool (const anArray: PDoubleArray;
  31. const X: WORD) : Double;
  32. begin
  33. Result := anArray^ [x];
  34. end;
  35.  
  36. function Frequency (const aValue: Double;
  37. const anArray: PDoubleArray;
  38. const aLength: WORD;
  39. const aPrecision: Double = 0.001) : Double;
  40. var
  41. I : WORD;
  42. Count : WORD;
  43. begin
  44. Result := 0;
  45. try
  46. Count := 0;
  47. for I := 0 to (aLength - 1) do
  48. begin
  49. if ( abs (aValue - anArray^ [I]) <= aPrecision ) then
  50. Inc (Count);
  51. Continue;
  52. end;
  53. Result := Count / aLength;
  54. except on E: Exception do
  55. raise Exception.CreateFmt ('%s : %s',[ERR_FREQUENCY,E.Message]);
  56. end;
  57. end;
  58.  
  59. function Entropy (const anArray: PDoubleArray;
  60. const aLength: WORD;
  61. const aPrecision: Double = 0.001) : Double;
  62. var
  63. I : WORD;
  64. Frq : Double;
  65. begin
  66. Result := 0;
  67. try
  68. for I := 0 to (aLength - 1) do
  69. begin
  70. Frq := Frequency (anArray^ [I],
  71. anArray,
  72. aLength,
  73. aPrecision);
  74. if ( Frq > 0 ) then
  75. Result := Result + Frq * log2 (1/Frq);
  76. Continue;
  77. end;
  78. except on E: Exception do
  79. raise Exception.CreateFmt ('%s : %s',[ERR_ENTROPY,E.Message]);
  80. end;
  81. end;
  82.  
  83. procedure FFT (var anArray: PDoubleArray;
  84. const aLength: WORD;
  85. const doInverse: Boolean = FALSE);
  86. var
  87. nn : WORD;
  88. ii : Integer;
  89. jj : Integer;
  90. n : Integer;
  91. m : Integer;
  92. mMax : Integer;
  93. i : Integer;
  94. j : Integer;
  95. iSign : Integer;
  96. iStep : Integer;
  97. wTemp : Double;
  98. wR : Double;
  99. wI : Double;
  100. wpR : Double;
  101. wpI : Double;
  102. tempR : Double;
  103. tempI : Double;
  104. theta : Double;
  105. begin
  106. try
  107. case doInverse of
  108. TRUE : iSign := 1;
  109. FALSE : iSign := -1;
  110. end;
  111. nn := aLength;
  112. n := 2*nn;
  113. j := 1;
  114. ii := 1;
  115. while ( ii <= nn ) do
  116. begin
  117. i := 2*ii - 1;
  118. if ( j > i ) then
  119. begin
  120. tempR := anArray^ [j-1];
  121. tempI := anArray^ [j];
  122. anArray^ [j-1] := anArray^ [i-1];
  123. anArray^ [j] := anArray^ [i];
  124. anArray^ [i-1] := tempR;
  125. anArray^ [i] := tempI;
  126. end;
  127. m := n div 2;
  128. while ( ( m >= 2) and ( j > m ) ) do
  129. begin
  130. j := j - m;
  131. m := m div 2;
  132. end;
  133. j := j + m;
  134. Inc (ii);
  135. end;
  136. mMax := 2;
  137. while ( n > mMax ) do
  138. begin
  139. iStep := 2*mMax;
  140. theta := 2 * Pi / (iSign*mMax);
  141. wpR := -2.0 * sqr ( sin (0.5*theta) );
  142. wpI := sin (theta);
  143. wR := 1.0;
  144. wI := 0.0;
  145. ii := 1;
  146. while ( ii <= mMax div 2 ) do
  147. begin
  148. m := 2*ii - 1;
  149. jj := 0;
  150. while ( jj <= (n-m) div iStep ) do
  151. begin
  152. i := m + jj*iStep;
  153. j := i + mMax;
  154. tempR := wR * anArray^ [j-1] - wI * anArray^ [j];
  155. tempI := wR * anArray^ [j] + wI * anArray^ [j-1];
  156. anArray^ [j-1] := anArray^ [i-1] - tempR;
  157. anArray^ [j] := anArray^ [i] - tempI;
  158. anArray^ [i-1] := anArray^ [i-1] + tempR;
  159. anArray^ [i] := anArray^ [i] + tempI;
  160. Inc (jj);
  161. end;
  162. wTemp := wR;
  163. wR := wR*wpR - wI*wpI + wR;
  164. wI := wI*wpR + wTemp*wpI + wI;
  165. Inc (ii);
  166. end;
  167. mMax := iStep;
  168. end;
  169. if doInverse then
  170. begin
  171. i := 1;
  172. while ( i <= 2*nn ) do
  173. begin
  174. anArray^ [i-1] := anArray^ [i-1] / nn;
  175. Inc (i);
  176. end;
  177. end;
  178. except on E: Exception do
  179. raise Exception.CreateFmt ('%s : %s',[ERR_FFT,E.Message]);
  180. end;
  181. end;
  182.  
  183. function DrawFunction (const aCanvas: TCanvas;
  184. const aRect: TRect;
  185. const aFunction: TDoubleFunction;
  186. const anArray: PDoubleArray;
  187. const aParams: TDoubleFunctionParams;
  188. const aBackGroundColor: TColor;
  189. const aBackColor: TColor;
  190. const aFrontColor: TColor;
  191. const drawAxes: Boolean = TRUE;
  192. const drawLines: Boolean = TRUE;
  193. const drawSubGraph: Boolean = FALSE) : Boolean;
  194. var
  195. x1, x2 : WORD; {ь значения функции }
  196.   x : WORD; { локальное значение аргумента }
  197.   y : Double; { локальное значение функц
  198. x : WORD; { { приращение аргумента }
  199.   x0, y0 : Integer; { начало координат }
  200.   Left : Integer;
  201.   Bottom : Integer;
  202.   Width : Integer;
  203.   Height : Integer;
  204. x0, y0 : Integer; { масштаб по оси x }
  205.   ScaleY : Double; { масштаб по оси y }
  206. begin
  207.   Result := FALSE;
  208.   try
  209.   x1 := aParams.X1;
  210.   x2 := aParams.X2;
  211.   dx := aPara
  212. ScaleY : Double; {vas do
  213.   begin
  214.   Brush.Color := aBackGroundColor;
  215.   Pen.Color := aFrontColor;
  216.   Rectangle (aRect);
  217.   if ( ( x2 < 0 ) or ( x1 < 0 ) ) then
  218.   Left := (aRect.Right - aRect.Left) div 2
  219.   else
  220.   Left := aRect.Left + 8;
  221.   Width := (aRect.Right - aRect.Left) - 18;
  222.   Bottom := (aRect.Bottom - aRect.Top) - 5;
  223.   Height := (aRect.Bottom - aRect.Top) - 22;
  224.   y1 := aFunction (anArray,x1);
  225.   y2 := aFunction (anArray,x1);
  226.   x := x1;
  227.   repeat
  228.   y := aFunction (anArray,x);
  229.   if ( y < y1 ) then
  230.   y1 := y;
  231.   if ( y > y2 ) then
  232.   y2 := y;
  233.   x := x + dx;
  234.   until ( x >= x2 );
  235.   if ( abs (y2-y1) < dx ) then Exit;
  236.   ScaleX := Width / abs (x2-x1);
  237.   ScaleY := Height / abs (y2-y1);
  238.   x0 := Left;
  239.   y0 := Bottom - abs ( round (y1*ScaleY) );
  240.   if drawAxes then
  241.   begin
  242.   MoveTo (Left,Bottom);
  243.   LineTo ( Left, Bottom-Height );
  244.   MoveTo (Left,y0);
  245.   LineTo ( x0+Width, y0 );
  246.   end;
  247.   Pen.Color := aBackColor;
  248.   x := x1;
  249.   repeat
  250.   y := aFunction (anArray,x);
  251.   if ( y > dx ) then
  252.   Pixels [ x0 + round (x*ScaleX), y0 - round (y*ScaleY)-1 ] := aFrontColor
  253.   else
  254.   Pixels [ x0 + round (x*ScaleX), y0 - round (y*ScaleY){+1} ] := aFrontColor;
  255.   MoveTo ( x0 + round (x*ScaleX), y0 - round (y*ScaleY) );
  256.   if not drawLines then
  257.   ] := aFrontColor;
  258. MoveTo ( x0 + round (x*ScaleX), y0 - round (y*ScaleY) );
  259. if not drawLines then
  260. x := x + dx
  261. else
  262. begin
  263. if drawSubGraph then
  264. begin
  265. LineTo ( x0 + round (x*ScaleX), y0 );
  266. x := x + dx;
  267. end
  268. else
  269. begin
  270. x := x + dx;
  271. y := aFunction (anArray,x);
  272. LineTo ( x0 + round (x*ScaleX), y0 - round (y*ScaleY) );
  273. end;
  274. end;
  275. until ( x >= x2 );
  276. Result := TRUE;
  277. end;
  278. except on E: Exception do
  279. raise Exception.CreateFmt ('%s : %s',[ERR_DRAW_FUNCTION,E.Message]);
  280. end;
  281. end;
  282.  
  283. //{ array [2n] + i*array [2n+1] }
  284. for I := 0 to ( 2 * f_PoolLength - 1 ) do
  285. begin
  286. f_Pool^ [I] := TCrypto.Random (f_PoolLength,0,f_RandomType);
  287. f_Pool^ [I
  288. for I := 0 to ( 2 * f_PoolLength - 1 ) do
  289. begin
  290. f_Pool^ [I] := TCrypto.Random (f_PoolLength,0,f_RandomType);
  291. f_Pool^ [I+1] := 0;
  292. end;
  293. Rect.Left := 4;
  294. Rect.Right := ptSpectr.Width - 4;
  295. Rect.Top := 4;
  296. Rect.Bottom := ptSpectr.Height - 4;
  297. Params.X1 := 1;
  298. Params.X2 := f_PoolLength;
  299. Params.dX := 2;
  300. DrawFunction (ptSpectr.Canvas,
  301. Rect,
  302. spectr_pool,
  303. f_Pool,
  304. Params,
  305. clBlack,
  306. clMaroon,
  307. clRed,
  308. TRUE,
  309. FALSE);
  310.  
  311. if f_EntropyPoolIndex >= f_EntropyPoolLength then
  312. f_EntropyPoolIndex := 0;
  313. Inc (f_EntropyPoolIndex);
  314. f_EntropyPool^ [f_EntropyPoolIndex] := Entropy ( f_Pool, 2*f_PoolLength );
  315. Rect.Left := 4;
  316. Rect.Right := ptEntropy.Width - 4;
  317. Rect.Top := 4;
  318. Rect.Bottom := ptEntropy.Height - 4;
  319. Params.X1 := 0;
  320. Params.X2 := f_EntropyPoolLength div 2;
  321. Params.dX := 1;
  322. DrawFunction (ptEntropy.Canvas,
  323. Rect,
  324. entropy_pool,
  325. f_EntropyPool,
  326. Params,
  327. clBlack,
  328. clMaroon,
  329. clRed,
  330. TRUE,
  331. TRUE,
  332. TRUE);
  333.  
  334. FFT (f_Pool,f_PoolLength);
  335. Rect.Left := 4;
  336. Rect.Right := ptFFT.Width - 4;
  337. Rect.Top := 4;
  338. Rect.Bottom := ptFFT.Height - 4;
  339. Params.X1 := 2;
  340. Params.X2 := f_PoolLength;
  341. Params.dX := 1;
  342. DrawFunction (ptFFT.Canvas,
  343. Rect,
  344. fft_pool,
  345. f_Pool,
  346. Params,
  347. clBlack,
  348. clMaroon,
  349. clRed,
  350. TRUE,
  351. TRUE,
  352. TRUE);
  353. except on E: Exception do
  354. raise

Ссылка на данный код:

На главную страницу сервиса обмена кодом »