Luxophia

1
unit LUX.DiscreteTrans.D1;
2

3
interface //#################################################################### ■
4

5
uses LUX;
6

7
type //$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【型】
8

9
     TDiscreteTrans = class;
10

11
     //$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【レコード】
12

13
     //$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【クラス】
14

15
     //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% TDiscreteTrans
16

17
     TDiscreteTrans = class
18
     private
19
       _NC :Integer;
20
       ///// initializing routines
21
       procedure makeipt;
22
       procedure makewt;
23
       procedure makect( var c:array of Double );
24
       /////
25
       procedure bitrv2;
26
       procedure bitrv2conj;
27

28
       procedure bitrv216;
29
       procedure bitrv216neg;
30

31
       procedure bitrv208;
32
       procedure bitrv208neg;
33

34
       procedure cftf1st( var w:array of Double );
35
       procedure cftb1st( var w:array of Double );
36

37
       procedure cftmdl1( const n:Integer; var a:array of Double; var w:array of Double );
38
       procedure cftmdl2( const n:Integer; var a:array of Double; var w:array of Double );
39

40
       function cfttree( const n,j,k:Integer ) :Integer;
41

42
       procedure cftf161( var a:array of Double; var w:array of Double );
43
       procedure cftf162( var a:array of Double; var w:array of Double );
44

45
       procedure cftf081( var a:array of Double; var w:array of Double );
46
       procedure cftf082( var a:array of Double; var w:array of Double );
47

48
       procedure cftleaf( const n,isplt:Integer; var a:array of Double );
49

50
       procedure cftrec4;
51

52
       procedure cftfx41;
53

54
       procedure cftf040;
55
       procedure cftb040;
56

57
       procedure cftx020;
58
     protected
59
       _TEMP  :array of Double;
60
       _TempN :Integer;
61
       _IP    :array of Integer;
62
       _W     :array of Double;
63
       _NW    :Integer;
64
       _NormW :Double;
65
       ///// フィールド
66
       _Count :Integer;  upCount :Boolean;
67
       ///// アクセス
68
       procedure SetCount( const Count_:Integer ); virtual;
69
       ///// メソッド
70
       procedure cftfsub;
71
       procedure cftbsub;
72
       procedure rftfsub( var c:array of Double );
73
       procedure rftbsub( var c:array of Double );
74

75
       procedure Normalize;
76

77
       procedure MakeTableW;
78
       procedure MakeTableC;
79
     public
80
       constructor Create;
81
       destructor Destroy; override;
82
       ///// プロパティ
83
       property Count :Integer read _Count write SetCount;
84
       ///// メソッド
85
       procedure TransWF; virtual; abstract;
86
       procedure TransFW; virtual; abstract;
87
     end;
88

89
//const //$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【定数】
90

91
//var //$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【変数】
92

93
//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【ルーチン】
94

95
implementation //################################################################################### ■
96

97
uses Math;
98

99
//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【レコード】
100

101
//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【クラス】
102

103
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% TDiscreteTrans
104

105
//&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& private
106

107
////////////////////////////////////////////////////////////////////////////// initializing routines
108

109
procedure TDiscreteTrans.makeipt;
110
var
111
   j, l, m, m2, p, q :Integer;
112
begin
113
     _IP[ 2 ] :=  0;
114
     _IP[ 3 ] := 16;
115

116
     m := 2;
117
     l := _NW;
118
     while l > 32 do
119
     begin
120
          m2 := m  shl 1;
121
          q  := m2 shl 3;
122

123
          for j := m to m2 - 1 do
124
          begin
125
               p := _IP[ j ] shl 2;
126

127
               _IP[ m  + j ] := p;
128
               _IP[ m2 + j ] := p + q;
129
          end;
130

131
          m := m2;
132
          l := l shr 2
133
     end
134
end;
135

136
procedure TDiscreteTrans.makewt;
137
var
138
   j, nwh, nw0, nw1 :Integer;
139
   delta, wn4r, wk1r, wk1i, wk3r, wk3i :Double;
140
begin
141
     _IP[ 0 ] := _NW;
142
     _IP[ 1 ] := 1;
143

144
     if _NW > 2 then
145
     begin
146
          nwh   := _NW shr 1;
147
          delta := ArcTan( 1 ) / nwh;
148
          wn4r  := Cos( delta * nwh );
149

150
          _W[ 0 ] := 1;
151
          _W[ 1 ] := wn4r;
152

153
          if nwh = 4 then
154
          begin
155
               _W[ 2 ] := Cos( delta * 2 );
156
               _W[ 3 ] := Sin( delta * 2 );
157
          end
158
          else
159
          if nwh > 4 then
160
          begin
161
               makeipt;
162

163
               _W[ 2 ] := 0.5 / Cos( delta * 2 );
164
               _W[ 3 ] := 0.5 / Cos( delta * 6 );
165

166
               j := 4;
167
               while j < nwh do
168
               begin
169
                    _W[ j     ] := +Cos(     delta * j );
170
                    _W[ j + 1 ] := +Sin(     delta * j );
171
                    _W[ j + 2 ] := +Cos( 3 * delta * j );
172
                    _W[ j + 3 ] := -Sin( 3 * delta * j );
173

174
                    Inc( j, 4 )
175
               end
176
          end;
177

178
          nw0 := 0;
179
          while nwh > 2 do
180
          begin
181
               nw1 := nw0 + nwh;
182
               nwh := nwh shr 1;
183

184
               _W[ nw1     ] := 1;
185
               _W[ nw1 + 1 ] := wn4r;
186

187
               if nwh = 4 then
188
               begin
189
                    wk1r := _W[ nw0 + 4 ];
190
                    wk1i := _W[ nw0 + 5 ];
191

192
                    _W[ nw1 + 2 ] := wk1r;
193
                    _W[ nw1 + 3 ] := wk1i;
194
               end
195
               else
196
               if nwh > 4 then
197
               begin
198
                    wk1r := _W[ nw0 + 4 ];
199
                    wk3r := _W[ nw0 + 6 ];
200

201
                    _W[ nw1 + 2 ] := 0.5 / wk1r;
202
                    _W[ nw1 + 3 ] := 0.5 / wk3r;
203

204
                    j := 4;
205
                    while j < nwh do
206
                    begin
207
                         wk1r := _W[ nw0 + 2 * j     ];
208
                         wk1i := _W[ nw0 + 2 * j + 1 ];
209
                         wk3r := _W[ nw0 + 2 * j + 2 ];
210
                         wk3i := _W[ nw0 + 2 * j + 3 ];
211

212
                         _W[ nw1 + j     ] := wk1r;
213
                         _W[ nw1 + j + 1 ] := wk1i;
214
                         _W[ nw1 + j + 2 ] := wk3r;
215
                         _W[ nw1 + j + 3 ] := wk3i;
216

217
                         Inc( j, 4 )
218
                    end
219
               end;
220

221
               nw0 := nw1
222
          end
223
     end
224
end;
225

226
procedure TDiscreteTrans.makect( var c:array of Double );
227
var
228
   j, nch :Integer;
229
   delta :Double;
230
begin
231
     _IP[ 1 ] := _NC;
232

233
     if _NC > 1 then
234
     begin
235
          nch   := _NC shr 1;
236
          delta := ArcTan( 1 ) / nch;
237

238
          c[  0  ] := Cos( delta * nch );
239
          c[ nch ] := 0.5 * c[ 0 ];
240

241
          for j := 1 to nch - 1do
242
          begin
243
               c[       j ] := 0.5 * Cos( delta * j );
244
               c[ _NC - j ] := 0.5 * Sin( delta * j );
245
          end
246
     end
247
end;
248

249
////////////////////////////////////////////////////////////////////////////////////////////////////
250

251
procedure TDiscreteTrans.bitrv2;
252
var
253
   j, j1, k, k1, l, m, nh, nm :Integer;
254
   xr, xi, yr, yi :Double;
255
begin
256
     m := 1;
257
     l := _TempN shr 2;
258
     while l > 8 do
259
     begin
260
          m := m shl 1;
261
          l := l shr 2
262
     end;
263

264
     nh := _TempN shr 1;
265
     nm := 4 * m;
266
     if l = 8 then
267
     begin
268
          for k := 0 to m - 1 do
269
          begin
270
               for j := 0 to k - 1 do
271
               begin
272
                    j1 := 4 * j + 2 * _IP[ m + k ];
273
                    k1 := 4 * k + 2 * _IP[ m + j ];
274

275
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
276
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
277
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
278
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
279

280
                    Inc( j1,     nm );
281
                    Inc( k1, 2 * nm );
282
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
283
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
284
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
285
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
286

287
                    Inc( j1, nm );
288
                    Dec( k1, nm );
289
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
290
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
291
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
292
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
293

294
                    Inc( j1,     nm );
295
                    Inc( k1, 2 * nm );
296
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
297
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
298
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
299
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
300

301
                    Inc( j1, nh );
302
                    Inc( k1,  2 );
303
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
304
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
305
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
306
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
307

308
                    Dec( j1,     nm );
309
                    Dec( k1, 2 * nm );
310
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
311
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
312
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
313
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
314

315
                    Dec( j1, nm );
316
                    Inc( k1, nm );
317
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
318
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
319
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
320
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
321

322
                    Dec( j1,     nm );
323
                    Dec( k1, 2 * nm );
324
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
325
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
326
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
327
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
328

329
                    Inc( j1,  2 );
330
                    Inc( k1, nh );
331
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
332
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
333
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
334
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
335

336
                    Inc( j1,     nm );
337
                    Inc( k1, 2 * nm );
338
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
339
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
340
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
341
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
342

343
                    Inc( j1, nm );
344
                    Dec( k1, nm );
345
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
346
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
347
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
348
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
349

350
                    Inc( j1,     nm );
351
                    Inc( k1, 2 * nm );
352
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
353
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
354
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
355
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
356

357
                    Dec( j1, nh );
358
                    Dec( k1,  2 );
359
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
360
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
361
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
362
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
363

364
                    Dec( j1, nm );
365
                    Dec( k1, 2 * nm );
366
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
367
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
368
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
369
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
370

371
                    Dec( j1, nm );
372
                    Inc( k1, nm );
373
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
374
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
375
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
376
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
377

378
                    Dec( j1, nm );
379
                    Dec( k1, 2 * nm );
380
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
381
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
382
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
383
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
384
               end;
385

386
               k1 := 4 * k + 2 * _IP[ m + k ];
387
               j1 := k1 + 2;
388
               Inc( k1, nh );
389
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
390
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
391
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
392
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
393

394
               Inc( j1, nm );
395
               Inc( k1, 2 * nm );
396
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
397
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
398
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
399
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
400

401
               Inc( j1, nm );
402
               Dec( k1, nm );
403
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
404
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
405
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
406
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
407

408
               Dec( j1, 2 );
409
               Dec( k1, nh );
410
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
411
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
412
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
413
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
414

415
               Inc( j1, nh + 2 );
416
               Inc( k1, nh + 2 );
417
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
418
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
419
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
420
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
421

422
               Dec( j1, nh - nm );
423
               Inc( k1, 2 * nm - 2 );
424
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
425
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
426
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
427
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
428
          end
429
     end
430
     else
431
     begin
432
          for k := 0 to m - 1 do
433
          begin
434
               for j := 0 to k - 1 do
435
               begin
436
                    j1 := 4 * j + _IP[ m + k ];
437
                    k1 := 4 * k + _IP[ m + j ];
438
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
439
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
440
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
441
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
442

443
                    Inc( j1, nm );
444
                    Inc( k1, nm );
445
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
446
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
447
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
448
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
449

450
                    Inc( j1, nh );
451
                    Inc( k1, 2 );
452
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
453
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
454
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
455
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
456

457
                    Dec( j1, nm );
458
                    Dec( k1, nm );
459
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
460
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
461
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
462
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
463

464
                    Inc( j1, 2 );
465
                    Inc( k1, nh );
466
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
467
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
468
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
469
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
470

471
                    Inc( j1, nm );
472
                    Inc( k1, nm );
473
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
474
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
475
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
476
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
477

478
                    Dec( j1, nh );
479
                    Dec( k1, 2 );
480
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
481
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
482
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
483
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
484

485
                    Dec( j1, nm );
486
                    Dec( k1, nm );
487
                    xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
488
                    yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
489
                    _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
490
                    _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
491
               end;
492

493
               k1 := 4 * k + _IP[ m + k ];
494
               j1 := k1 + 2;
495
               Inc( k1, nh );
496
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
497
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
498
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
499
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
500

501
               Inc( j1, nm );
502
               Inc( k1, nm );
503
               xr := _TEMP[ j1 ]; xi := _TEMP[ j1 + 1 ];
504
               yr := _TEMP[ k1 ]; yi := _TEMP[ k1 + 1 ];
505
               _TEMP[ j1 ] := yr; _TEMP[ j1 + 1 ] := yi;
506
               _TEMP[ k1 ] := xr; _TEMP[ k1 + 1 ] := xi;
507
          end
508
     end
509
end;
510

511
procedure TDiscreteTrans.bitrv2conj;
512
var
513
   j, j1, k, k1, l, m, nh, nm :Integer;
514
   xr, xi, yr, yi :Double;
515
begin
516
     m := 1;
517
     l := _TempN shr 2;
518
     while l > 8 do
519
     begin
520
          m := m shl 1;
521
          l := l shr 2;
522
     end;
523

524
     nh := _TempN shr 1;
525
     nm := 4 * m;
526
     if l = 8 then
527
     begin
528
          for k := 0 to m - 1 do
529
          begin
530
               for j := 0 to k - 1 do
531
               begin
532
                    j1 := 4 * j + 2 * _IP[ m + k ];
533
                    k1 := 4 * k + 2 * _IP[ m + j ];
534

535
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
536
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
537
                    _TEMP[ j1 ] := yr;
538
                    _TEMP[ j1 + 1 ] := yi;
539
                    _TEMP[ k1 ] := xr;
540
                    _TEMP[ k1 + 1 ] := xi;
541

542
                    Inc( j1,     nm );
543
                    Inc( k1, 2 * nm );
544
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
545
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
546
                    _TEMP[ j1     ] := yr;
547
                    _TEMP[ j1 + 1 ] := yi;
548
                    _TEMP[ k1     ] := xr;
549
                    _TEMP[ k1 + 1 ] := xi;
550

551
                    Inc( j1, nm );
552
                    Dec( k1, nm );
553
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
554
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
555
                    _TEMP[ j1     ] := yr;
556
                    _TEMP[ j1 + 1 ] := yi;
557
                    _TEMP[ k1     ] := xr;
558
                    _TEMP[ k1 + 1 ] := xi;
559

560
                    Inc( j1,     nm );
561
                    Inc( k1, 2 * nm );
562
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
563
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
564
                    _TEMP[ j1     ] := yr;
565
                    _TEMP[ j1 + 1 ] := yi;
566
                    _TEMP[ k1     ] := xr;
567
                    _TEMP[ k1 + 1 ] := xi;
568

569
                    Inc( j1, nh );
570
                    Inc( k1,  2 );
571
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
572
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
573
                    _TEMP[ j1     ] := yr;
574
                    _TEMP[ j1 + 1 ] := yi;
575
                    _TEMP[ k1     ] := xr;
576
                    _TEMP[ k1 + 1 ] := xi;
577

578
                    Dec( j1,     nm );
579
                    Dec( k1, 2 * nm );
580
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
581
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
582
                    _TEMP[ j1     ] := yr;
583
                    _TEMP[ j1 + 1 ] := yi;
584
                    _TEMP[ k1     ] := xr;
585
                    _TEMP[ k1 + 1 ] := xi;
586

587
                    Dec( j1, nm );
588
                    Inc( k1, nm );
589
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
590
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
591
                    _TEMP[ j1     ] := yr;
592
                    _TEMP[ j1 + 1 ] := yi;
593
                    _TEMP[ k1     ] := xr;
594
                    _TEMP[ k1 + 1 ] := xi;
595

596
                    Dec( j1,     nm );
597
                    Dec( k1, 2 * nm );
598
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
599
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
600
                    _TEMP[ j1     ] := yr;
601
                    _TEMP[ j1 + 1 ] := yi;
602
                    _TEMP[ k1     ] := xr;
603
                    _TEMP[ k1 + 1 ] := xi;
604

605
                    Inc( j1,  2 );
606
                    Inc( k1, nh );
607
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
608
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
609
                    _TEMP[ j1     ] := yr;
610
                    _TEMP[ j1 + 1 ] := yi;
611
                    _TEMP[ k1     ] := xr;
612
                    _TEMP[ k1 + 1 ] := xi;
613

614
                    Inc( j1,     nm );
615
                    Inc( k1, 2 * nm );
616
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
617
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
618
                    _TEMP[ j1     ] := yr;
619
                    _TEMP[ j1 + 1 ] := yi;
620
                    _TEMP[ k1     ] := xr;
621
                    _TEMP[ k1 + 1 ] := xi;
622

623
                    Inc( j1, nm );
624
                    Dec( k1, nm );
625
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
626
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
627
                    _TEMP[ j1     ] := yr;
628
                    _TEMP[ j1 + 1 ] := yi;
629
                    _TEMP[ k1     ] := xr;
630
                    _TEMP[ k1 + 1 ] := xi;
631

632
                    Inc( j1,     nm );
633
                    Inc( k1, 2 * nm );
634
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
635
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
636
                    _TEMP[ j1     ] := yr;
637
                    _TEMP[ j1 + 1 ] := yi;
638
                    _TEMP[ k1     ] := xr;
639
                    _TEMP[ k1 + 1 ] := xi;
640

641
                    Dec( j1, nh );
642
                    Dec( k1,  2 );
643
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
644
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
645
                    _TEMP[ j1     ] := yr;
646
                    _TEMP[ j1 + 1 ] := yi;
647
                    _TEMP[ k1     ] := xr;
648
                    _TEMP[ k1 + 1 ] := xi;
649

650
                    Dec( j1,     nm );
651
                    Dec( k1, 2 * nm );
652
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
653
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
654
                    _TEMP[ j1     ] := yr;
655
                    _TEMP[ j1 + 1 ] := yi;
656
                    _TEMP[ k1     ] := xr;
657
                    _TEMP[ k1 + 1 ] := xi;
658

659
                    Dec( j1, nm );
660
                    Inc( k1, nm );
661
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
662
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
663
                    _TEMP[ j1     ] := yr;
664
                    _TEMP[ j1 + 1 ] := yi;
665
                    _TEMP[ k1     ] := xr;
666
                    _TEMP[ k1 + 1 ] := xi;
667

668
                    Dec( j1,     nm );
669
                    Dec( k1, 2 * nm );
670
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
671
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
672
                    _TEMP[ j1     ] := yr;
673
                    _TEMP[ j1 + 1 ] := yi;
674
                    _TEMP[ k1     ] := xr;
675
                    _TEMP[ k1 + 1 ] := xi;
676
               end;
677

678
               k1 := 4 * k + 2 * _IP[ m + k ];
679

680
               j1 := k1 + 2;
681
               Inc( k1, nh );
682
               _TEMP[ j1 - 1 ] := -_TEMP[ j1 - 1 ];
683
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
684
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
685
               _TEMP[ j1     ] := yr;
686
               _TEMP[ j1 + 1 ] := yi;
687
               _TEMP[ k1     ] := xr;
688
               _TEMP[ k1 + 1 ] := xi;
689
               _TEMP[ k1 + 3 ] := -_TEMP[ k1 + 3 ];
690

691
               Inc( j1, nm );
692
               Inc( k1, 2 * nm );
693
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
694
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
695
               _TEMP[ j1     ] := yr;
696
               _TEMP[ j1 + 1 ] := yi;
697
               _TEMP[ k1     ] := xr;
698
               _TEMP[ k1 + 1 ] := xi;
699

700
               Inc( j1, nm );
701
               Dec( k1, nm );
702
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
703
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
704
               _TEMP[ j1     ] := yr;
705
               _TEMP[ j1 + 1 ] := yi;
706
               _TEMP[ k1     ] := xr;
707
               _TEMP[ k1 + 1 ] := xi;
708

709
               Dec( j1, 2 );
710
               Dec( k1, nh );
711
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
712
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
713
               _TEMP[ j1     ] := yr;
714
               _TEMP[ j1 + 1 ] := yi;
715
               _TEMP[ k1     ] := xr;
716
               _TEMP[ k1 + 1 ] := xi;
717

718
               Inc( j1, nh + 2 );
719
               Inc( k1, nh + 2 );
720
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
721
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
722
               _TEMP[ j1     ] := yr;
723
               _TEMP[ j1 + 1 ] := yi;
724
               _TEMP[ k1    ] := xr;
725
               _TEMP[ k1 + 1 ] := xi;
726

727
               Dec( j1, nh - nm );
728
               Inc( k1, 2 * nm - 2 );
729
               _TEMP[ j1 - 1 ] := -_TEMP[ j1 - 1 ];
730
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
731
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
732
               _TEMP[ j1     ] := yr;
733
               _TEMP[ j1 + 1 ] := yi;
734
               _TEMP[ k1     ] := xr;
735
               _TEMP[ k1 + 1 ] := xi;
736
               _TEMP[ k1 + 3 ] := -_TEMP[ k1 + 3 ];
737
          end
738
     end
739
     else
740
     begin
741
          for k := 0 to m - 1 do
742
          begin
743
               for j := 0 to k - 1 do
744
               begin
745
                    j1 := 4 * j + _IP[ m + k ];
746
                    k1 := 4 * k + _IP[ m + j ];
747

748
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
749
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
750
                    _TEMP[ j1     ] := yr;
751
                    _TEMP[ j1 + 1 ] := yi;
752
                    _TEMP[ k1     ] := xr;
753
                    _TEMP[ k1 + 1 ] := xi;
754

755
                    Inc( j1, nm );
756
                    Inc( k1, nm );
757
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
758
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
759
                    _TEMP[ j1     ] := yr;
760
                    _TEMP[ j1 + 1 ] := yi;
761
                    _TEMP[ k1     ] := xr;
762
                    _TEMP[ k1 + 1 ] := xi;
763

764
                    Inc( j1, nh );
765
                    Inc( k1,  2 );
766
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
767
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
768
                    _TEMP[ j1     ] := yr;
769
                    _TEMP[ j1 + 1 ] := yi;
770
                    _TEMP[ k1     ] := xr;
771
                    _TEMP[ k1 + 1 ] := xi;
772

773
                    Dec( j1, nm );
774
                    Dec( k1, nm );
775
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
776
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
777
                    _TEMP[ j1     ] := yr;
778
                    _TEMP[ j1 + 1 ] := yi;
779
                    _TEMP[ k1     ] := xr;
780
                    _TEMP[ k1 + 1 ] := xi;
781

782
                    Inc( j1,  2 );
783
                    Inc( k1, nh );
784
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
785
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
786
                    _TEMP[ j1     ] := yr;
787
                    _TEMP[ j1 + 1 ] := yi;
788
                    _TEMP[ k1     ] := xr;
789
                    _TEMP[ k1 + 1 ] := xi;
790

791
                    Inc( j1, nm );
792
                    Inc( k1, nm );
793
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
794
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
795
                    _TEMP[ j1     ] := yr;
796
                    _TEMP[ j1 + 1 ] := yi;
797
                    _TEMP[ k1     ] := xr;
798
                    _TEMP[ k1 + 1 ] := xi;
799

800
                    Dec( j1, nh );
801
                    Dec( k1,  2 );
802
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
803
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
804
                    _TEMP[ j1     ] := yr;
805
                    _TEMP[ j1 + 1 ] := yi;
806
                    _TEMP[ k1     ] := xr;
807
                    _TEMP[ k1 + 1 ] := xi;
808

809
                    Dec( j1, nm );
810
                    Dec( k1, nm );
811
                    xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
812
                    yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
813
                    _TEMP[ j1     ] := yr;
814
                    _TEMP[ j1 + 1 ] := yi;
815
                    _TEMP[ k1     ] := xr;
816
                    _TEMP[ k1 + 1 ] := xi;
817
               end;
818

819
               k1 := 4 * k + _IP[ m + k ];
820
               j1 := k1 + 2;
821
               Inc( k1, nh );
822
               _TEMP[ j1 - 1 ] := -_TEMP[ j1 - 1 ];
823
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
824
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
825
               _TEMP[ j1     ] := yr;
826
               _TEMP[ j1 + 1 ] := yi;
827
               _TEMP[ k1     ] := xr;
828
               _TEMP[ k1 + 1 ] := xi;
829
               _TEMP[ k1 + 3 ] := -_TEMP[ k1 + 3 ];
830

831
               Inc( j1, nm );
832
               Inc( k1, nm );
833
               _TEMP[ j1 - 1 ] := -_TEMP[ j1 - 1 ];
834
               xr := +_TEMP[ j1     ]; xi := -_TEMP[ j1 + 1 ];
835
               yr := +_TEMP[ k1     ]; yi := -_TEMP[ k1 + 1 ];
836
               _TEMP[ j1     ] := yr;
837
               _TEMP[ j1 + 1 ] := yi;
838
               _TEMP[ k1     ] := xr;
839
               _TEMP[ k1 + 1 ] := xi;
840
               _TEMP[ k1 + 3 ] := -_TEMP[ k1 + 3 ];
841
          end
842
     end
843
end;
844

845
////////////////////////////////////////////////////////////////////////////////////////////////////
846

847
procedure TDiscreteTrans.bitrv216;
848
var
849
   x01r, x01i,
850
   x02r, x02i,
851
   x03r, x03i,
852
   x04r, x04i,
853
   x05r, x05i,
854
   {--}  {--}
855
   x07r, x07i,
856
   x08r, x08i,
857
   {--}  {--}
858
   x10r, x10i,
859
   x11r, x11i,
860
   x12r, x12i,
861
   x13r, x13i,
862
   x14r, x14i
863
   {--}  {--} :Double;
864
begin
865
    x01r := _TEMP[ 02 ];  x01i := _TEMP[ 03 ];
866
    x02r := _TEMP[ 04 ];  x02i := _TEMP[ 05 ];
867
    x03r := _TEMP[ 06 ];  x03i := _TEMP[ 07 ];
868
    x04r := _TEMP[ 08 ];  x04i := _TEMP[ 09 ];
869
    x05r := _TEMP[ 10 ];  x05i := _TEMP[ 11 ];
870
    {-----------------}   {-----------------}
871
    x07r := _TEMP[ 14 ];  x07i := _TEMP[ 15 ];
872
    x08r := _TEMP[ 16 ];  x08i := _TEMP[ 17 ];
873
    {-----------------}   {-----------------}
874
    x10r := _TEMP[ 20 ];  x10i := _TEMP[ 21 ];
875
    x11r := _TEMP[ 22 ];  x11i := _TEMP[ 23 ];
876
    x12r := _TEMP[ 24 ];  x12i := _TEMP[ 25 ];
877
    x13r := _TEMP[ 26 ];  x13i := _TEMP[ 27 ];
878
    x14r := _TEMP[ 28 ];  x14i := _TEMP[ 29 ];
879
    {-----------------}   {-----------------}
880

881
    _TEMP[ 02 ] := x08r;  _TEMP[ 03 ] := x08i;
882
    _TEMP[ 04 ] := x04r;  _TEMP[ 05 ] := x04i;
883
    _TEMP[ 06 ] := x12r;  _TEMP[ 07 ] := x12i;
884
    _TEMP[ 08 ] := x02r;  _TEMP[ 09 ] := x02i;
885
    _TEMP[ 10 ] := x10r;  _TEMP[ 11 ] := x10i;
886
    {-----------------}   {-----------------}
887
    _TEMP[ 14 ] := x14r;  _TEMP[ 15 ] := x14i;
888
    _TEMP[ 16 ] := x01r;  _TEMP[ 17 ] := x01i;
889
    {-----------------}   {-----------------}
890
    _TEMP[ 20 ] := x05r;  _TEMP[ 21 ] := x05i;
891
    _TEMP[ 22 ] := x13r;  _TEMP[ 23 ] := x13i;
892
    _TEMP[ 24 ] := x03r;  _TEMP[ 25 ] := x03i;
893
    _TEMP[ 26 ] := x11r;  _TEMP[ 27 ] := x11i;
894
    _TEMP[ 28 ] := x07r;  _TEMP[ 29 ] := x07i;
895
    {-----------------}   {-----------------}
896
end;
897

898

899
procedure TDiscreteTrans.bitrv216neg;
900
var
901
   x01r, x01i,
902
   x02r, x02i,
903
   x03r, x03i,
904
   x04r, x04i,
905
   x05r, x05i,
906
   x06r, x06i,
907
   x07r, x07i,
908
   x08r, x08i,
909
   x09r, x09i,
910
   x10r, x10i,
911
   x11r, x11i,
912
   x12r, x12i,
913
   x13r, x13i,
914
   x14r, x14i,
915
   x15r, x15i :Double;
916
begin
917
     x01r := _TEMP[ 02 ];  x01i := _TEMP[ 03 ];
918
     x02r := _TEMP[ 04 ];  x02i := _TEMP[ 05 ];
919
     x03r := _TEMP[ 06 ];  x03i := _TEMP[ 07 ];
920
     x04r := _TEMP[ 08 ];  x04i := _TEMP[ 09 ];
921
     x05r := _TEMP[ 10 ];  x05i := _TEMP[ 11 ];
922
     x06r := _TEMP[ 12 ];  x06i := _TEMP[ 13 ];
923
     x07r := _TEMP[ 14 ];  x07i := _TEMP[ 15 ];
924
     x08r := _TEMP[ 16 ];  x08i := _TEMP[ 17 ];
925
     x09r := _TEMP[ 18 ];  x09i := _TEMP[ 19 ];
926
     x10r := _TEMP[ 20 ];  x10i := _TEMP[ 21 ];
927
     x11r := _TEMP[ 22 ];  x11i := _TEMP[ 23 ];
928
     x12r := _TEMP[ 24 ];  x12i := _TEMP[ 25 ];
929
     x13r := _TEMP[ 26 ];  x13i := _TEMP[ 27 ];
930
     x14r := _TEMP[ 28 ];  x14i := _TEMP[ 29 ];
931
     x15r := _TEMP[ 30 ];  x15i := _TEMP[ 31 ];
932

933
     _TEMP[ 02 ] := x15r;  _TEMP[ 03 ] := x15i;
934
     _TEMP[ 04 ] := x07r;  _TEMP[ 05 ] := x07i;
935
     _TEMP[ 06 ] := x11r;  _TEMP[ 07 ] := x11i;
936
     _TEMP[ 08 ] := x03r;  _TEMP[ 09 ] := x03i;
937
     _TEMP[ 10 ] := x13r;  _TEMP[ 11 ] := x13i;
938
     _TEMP[ 12 ] := x05r;  _TEMP[ 13 ] := x05i;
939
     _TEMP[ 14 ] := x09r;  _TEMP[ 15 ] := x09i;
940
     _TEMP[ 16 ] := x01r;  _TEMP[ 17 ] := x01i;
941
     _TEMP[ 18 ] := x14r;  _TEMP[ 19 ] := x14i;
942
     _TEMP[ 20 ] := x06r;  _TEMP[ 21 ] := x06i;
943
     _TEMP[ 22 ] := x10r;  _TEMP[ 23 ] := x10i;
944
     _TEMP[ 24 ] := x02r;  _TEMP[ 25 ] := x02i;
945
     _TEMP[ 26 ] := x12r;  _TEMP[ 27 ] := x12i;
946
     _TEMP[ 28 ] := x04r;  _TEMP[ 29 ] := x04i;
947
     _TEMP[ 30 ] := x08r;  _TEMP[ 31 ] := x08i;
948
end;
949

950
////////////////////////////////////////////////////////////////////////////////////////////////////
951

952
procedure TDiscreteTrans.bitrv208;
953
var
954
   x1r, x1i,
955
   {-}  {-}
956
   x3r, x3i,
957
   x4r, x4i,
958
   {-}  {-}
959
   x6r, x6i :Double;
960
begin
961
     x1r := _TEMP[ 02 ];  x1i := _TEMP[ 03 ];
962
     {----------------}   {----------------}
963
     x3r := _TEMP[ 06 ];  x3i := _TEMP[ 07 ];
964
     x4r := _TEMP[ 08 ];  x4i := _TEMP[ 09 ];
965
     {----------------}   {----------------}
966
     x6r := _TEMP[ 12 ];  x6i := _TEMP[ 13 ];
967

968
     _TEMP[ 02 ] := x4r;  _TEMP[ 03 ] := x4i;
969
     {----------------}   {----------------}
970
     _TEMP[ 06 ] := x6r;  _TEMP[ 07 ] := x6i;
971
     _TEMP[ 08 ] := x1r;  _TEMP[ 09 ] := x1i;
972
     {----------------}   {----------------}
973
     _TEMP[ 12 ] := x3r;  _TEMP[ 13 ] := x3i;
974
end;
975

976
procedure TDiscreteTrans.bitrv208neg;
977
var
978
   x1r, x1i,
979
   x2r, x2i,
980
   x3r, x3i,
981
   x4r, x4i,
982
   x5r, x5i,
983
   x6r, x6i,
984
   x7r, x7i :Double;
985
begin
986
     x1r := _TEMP[ 02 ];  x1i := _TEMP[ 03 ];
987
     x2r := _TEMP[ 04 ];  x2i := _TEMP[ 05 ];
988
     x3r := _TEMP[ 06 ];  x3i := _TEMP[ 07 ];
989
     x4r := _TEMP[ 08 ];  x4i := _TEMP[ 09 ];
990
     x5r := _TEMP[ 10 ];  x5i := _TEMP[ 11 ];
991
     x6r := _TEMP[ 12 ];  x6i := _TEMP[ 13 ];
992
     x7r := _TEMP[ 14 ];  x7i := _TEMP[ 15 ];
993

994
     _TEMP[ 02 ] := x7r;  _TEMP[ 03 ] := x7i;
995
     _TEMP[ 04 ] := x3r;  _TEMP[ 05 ] := x3i;
996
     _TEMP[ 06 ] := x5r;  _TEMP[ 07 ] := x5i;
997
     _TEMP[ 08 ] := x1r;  _TEMP[ 09 ] := x1i;
998
     _TEMP[ 10 ] := x6r;  _TEMP[ 11 ] := x6i;
999
     _TEMP[ 12 ] := x2r;  _TEMP[ 13 ] := x2i;
1000
     _TEMP[ 14 ] := x4r;  _TEMP[ 15 ] := x4i;
1001
end;
1002

1003
////////////////////////////////////////////////////////////////////////////////////////////////////
1004

1005
procedure TDiscreteTrans.cftf1st( var w:array of Double );
1006
var
1007
   j, j0, j1, j2, j3, k, m, mh :Integer;
1008
   wn4r, csc1, csc3, wk1r, wk1i, wk3r, wk3i,
1009
   wd1r, wd1i, wd3r, wd3i,
1010
   x0r, x0i, x1r, x1i, x2r, x2i, x3r, x3i,
1011
   y0r, y0i, y1r, y1i, y2r, y2i, y3r, y3i :Double;
1012
begin
1013
     mh := _TempN shr 3;
1014
     m := 2 * mh;
1015

1016
     j1 :=      m;
1017
     j2 := j1 + m;
1018
     j3 := j2 + m;
1019

1020
     x0r := _TEMP[ 0  ] + _TEMP[ j2 ]; x0i := _TEMP[ 1      ] + _TEMP[ j2 + 1 ];
1021
     x1r := _TEMP[ 0  ] - _TEMP[ j2 ]; x1i := _TEMP[ 1      ] - _TEMP[ j2 + 1 ];
1022
     x2r := _TEMP[ j1 ] + _TEMP[ j3 ]; x2i := _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1023
     x3r := _TEMP[ j1 ] - _TEMP[ j3 ]; x3i := _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1024
     _TEMP[ 0 ] := x0r + x2r;
1025
     _TEMP[ 1 ] := x0i + x2i;
1026
     _TEMP[ j1 ] := x0r - x2r;
1027
     _TEMP[ j1 + 1 ] := x0i - x2i;
1028
     _TEMP[ j2 ] := x1r - x3i;
1029
     _TEMP[ j2 + 1 ] := x1i + x3r;
1030
     _TEMP[ j3 ] := x1r + x3i;
1031
     _TEMP[ j3 + 1 ] := x1i - x3r;
1032

1033
     wn4r := w[ 1 ];
1034
     csc1 := w[ 2 ];
1035
     csc3 := w[ 3 ];
1036
     wd1r := 1;
1037
     wd1i := 0;
1038
     wd3r := 1;
1039
     wd3i := 0;
1040
     k := 0;
1041
     j := 2;
1042
     while j < mh - 2 do
1043
     begin
1044
          Inc( k, 4 );
1045

1046
          wk1r := csc1 * ( wd1r + w[ k     ] ); wk1i := csc1 * ( wd1i + w[ k + 1 ] );
1047
          wk3r := csc3 * ( wd3r + w[ k + 2 ] ); wk3i := csc3 * ( wd3i + w[ k + 3 ] );
1048

1049
          wd1r := w[ k     ]; wd1i := w[ k + 1 ];
1050
          wd3r := w[ k + 2 ]; wd3i := w[ k + 3 ];
1051

1052
          j1 := j  + m;
1053
          j2 := j1 + m;
1054
          j3 := j2 + m;
1055

1056
          x0r := _TEMP[ j      ] + _TEMP[ j2     ]; x0i := _TEMP[ j  + 1 ] + _TEMP[ j2 + 1 ];
1057
          x1r := _TEMP[ j      ] - _TEMP[ j2     ]; x1i := _TEMP[ j  + 1 ] - _TEMP[ j2 + 1 ];
1058
          y0r := _TEMP[ j  + 2 ] + _TEMP[ j2 + 2 ]; y0i := _TEMP[ j  + 3 ] + _TEMP[ j2 + 3 ];
1059
          y1r := _TEMP[ j  + 2 ] - _TEMP[ j2 + 2 ]; y1i := _TEMP[ j  + 3 ] - _TEMP[ j2 + 3 ];
1060
          x2r := _TEMP[ j1     ] + _TEMP[ j3     ]; x2i := _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1061
          x3r := _TEMP[ j1     ] - _TEMP[ j3     ]; x3i := _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1062
          y2r := _TEMP[ j1 + 2 ] + _TEMP[ j3 + 2 ]; y2i := _TEMP[ j1 + 3 ] + _TEMP[ j3 + 3 ];
1063
          y3r := _TEMP[ j1 + 2 ] - _TEMP[ j3 + 2 ]; y3i := _TEMP[ j1 + 3 ] - _TEMP[ j3 + 3 ];
1064
          _TEMP[ j ] := x0r + x2r;
1065
          _TEMP[ j + 1 ] := x0i + x2i;
1066
          _TEMP[ j + 2 ] := y0r + y2r;
1067
          _TEMP[ j + 3 ] := y0i + y2i;
1068
          _TEMP[ j1 ] := x0r - x2r;
1069
          _TEMP[ j1 + 1 ] := x0i - x2i;
1070
          _TEMP[ j1 + 2 ] := y0r - y2r;
1071
          _TEMP[ j1 + 3 ] := y0i - y2i;
1072

1073
          x0r := x1r - x3i; x0i := x1i + x3r;
1074
          _TEMP[ j2     ] := wk1r * x0r - wk1i * x0i;
1075
          _TEMP[ j2 + 1 ] := wk1r * x0i + wk1i * x0r;
1076

1077
          x0r := y1r - y3i; x0i := y1i + y3r;
1078
          _TEMP[ j2 + 2 ] := wd1r * x0r - wd1i * x0i;
1079
          _TEMP[ j2 + 3 ] := wd1r * x0i + wd1i * x0r;
1080

1081
          x0r := x1r + x3i; x0i := x1i - x3r;
1082
          _TEMP[ j3 ] := wk3r * x0r + wk3i * x0i;
1083
          _TEMP[ j3 + 1 ] := wk3r * x0i - wk3i * x0r;
1084

1085
          x0r := y1r + y3i; x0i := y1i - y3r;
1086
          _TEMP[ j3 + 2 ] := wd3r * x0r + wd3i * x0i;
1087
          _TEMP[ j3 + 3 ] := wd3r * x0i - wd3i * x0r;
1088

1089
          j0 :=  m - j;
1090
          j1 := j0 + m;
1091
          j2 := j1 + m;
1092
          j3 := j2 + m;
1093

1094
          x0r := _TEMP[ j0     ] + _TEMP[ j2     ]; x0i := _TEMP[ j0 + 1 ] + _TEMP[ j2 + 1 ];
1095
          x1r := _TEMP[ j0     ] - _TEMP[ j2     ]; x1i := _TEMP[ j0 + 1 ] - _TEMP[ j2 + 1 ];
1096
          y0r := _TEMP[ j0 - 2 ] + _TEMP[ j2 - 2 ]; y0i := _TEMP[ j0 - 1 ] + _TEMP[ j2 - 1 ];
1097
          y1r := _TEMP[ j0 - 2 ] - _TEMP[ j2 - 2 ]; y1i := _TEMP[ j0 - 1 ] - _TEMP[ j2 - 1 ];
1098
          x2r := _TEMP[ j1     ] + _TEMP[ j3     ]; x2i := _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1099
          x3r := _TEMP[ j1     ] - _TEMP[ j3     ]; x3i := _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1100
          y2r := _TEMP[ j1 - 2 ] + _TEMP[ j3 - 2 ]; y2i := _TEMP[ j1 - 1 ] + _TEMP[ j3 - 1 ];
1101
          y3r := _TEMP[ j1 - 2 ] - _TEMP[ j3 - 2 ]; y3i := _TEMP[ j1 - 1 ] - _TEMP[ j3 - 1 ];
1102
          _TEMP[ j0     ] := x0r + x2r;
1103
          _TEMP[ j0 + 1 ] := x0i + x2i;
1104
          _TEMP[ j0 - 2 ] := y0r + y2r;
1105
          _TEMP[ j0 - 1 ] := y0i + y2i;
1106
          _TEMP[ j1     ] := x0r - x2r;
1107
          _TEMP[ j1 + 1 ] := x0i - x2i;
1108
          _TEMP[ j1 - 2 ] := y0r - y2r;
1109
          _TEMP[ j1 - 1 ] := y0i - y2i;
1110

1111
          x0r := x1r - x3i; x0i := x1i + x3r;
1112
          _TEMP[ j2     ] := wk1i * x0r - wk1r * x0i;
1113
          _TEMP[ j2 + 1 ] := wk1i * x0i + wk1r * x0r;
1114

1115
          x0r := y1r - y3i; x0i := y1i + y3r;
1116
          _TEMP[ j2 - 2 ] := wd1i * x0r - wd1r * x0i;
1117
          _TEMP[ j2 - 1 ] := wd1i * x0i + wd1r * x0r;
1118

1119
          x0r := x1r + x3i; x0i := x1i - x3r;
1120
          _TEMP[ j3     ] := wk3i * x0r + wk3r * x0i;
1121
          _TEMP[ j3 + 1 ] := wk3i * x0i - wk3r * x0r;
1122

1123
          x0r := y1r + y3i; x0i := y1i - y3r;
1124
          _TEMP[ j3 - 2 ] := wd3i * x0r + wd3r * x0i;
1125
          _TEMP[ j3 - 1 ] := wd3i * x0i - wd3r * x0r;
1126

1127
          Inc( j, 4 )
1128
     end;
1129

1130
     wk1r := csc1 * ( wd1r + wn4r );
1131
     wk1i := csc1 * ( wd1i + wn4r );
1132
     wk3r := csc3 * ( wd3r - wn4r );
1133
     wk3i := csc3 * ( wd3i - wn4r );
1134

1135
     j0 := mh;
1136
     j1 := j0 + m;
1137
     j2 := j1 + m;
1138
     j3 := j2 + m;
1139

1140
     x0r := _TEMP[ j0 - 2 ] + _TEMP[ j2 - 2 ]; x0i := _TEMP[ j0 - 1 ] + _TEMP[ j2 - 1 ];
1141
     x1r := _TEMP[ j0 - 2 ] - _TEMP[ j2 - 2 ]; x1i := _TEMP[ j0 - 1 ] - _TEMP[ j2 - 1 ];
1142
     x2r := _TEMP[ j1 - 2 ] + _TEMP[ j3 - 2 ]; x2i := _TEMP[ j1 - 1 ] + _TEMP[ j3 - 1 ];
1143
     x3r := _TEMP[ j1 - 2 ] - _TEMP[ j3 - 2 ]; x3i := _TEMP[ j1 - 1 ] - _TEMP[ j3 - 1 ];
1144
     _TEMP[ j0 - 2 ] := x0r + x2r;
1145
     _TEMP[ j0 - 1 ] := x0i + x2i;
1146
     _TEMP[ j1 - 2 ] := x0r - x2r;
1147
     _TEMP[ j1 - 1 ] := x0i - x2i;
1148

1149
     x0r := x1r - x3i; x0i := x1i + x3r;
1150
     _TEMP[ j2 - 2 ] := wk1r * x0r - wk1i * x0i;
1151
     _TEMP[ j2 - 1 ] := wk1r * x0i + wk1i * x0r;
1152

1153
     x0r := x1r + x3i; x0i := x1i - x3r;
1154
     _TEMP[ j3 - 2 ] := wk3r * x0r + wk3i * x0i;
1155
     _TEMP[ j3 - 1 ] := wk3r * x0i - wk3i * x0r;
1156

1157
     x0r := _TEMP[ j0     ] + _TEMP[ j2     ]; x0i := _TEMP[ j0 + 1 ] + _TEMP[ j2 + 1 ];
1158
     x1r := _TEMP[ j0     ] - _TEMP[ j2     ]; x1i := _TEMP[ j0 + 1 ] - _TEMP[ j2 + 1 ];
1159
     x2r := _TEMP[ j1     ] + _TEMP[ j3     ]; x2i := _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1160
     x3r := _TEMP[ j1     ] - _TEMP[ j3     ]; x3i := _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1161
     _TEMP[ j0     ] := x0r + x2r;
1162
     _TEMP[ j0 + 1 ] := x0i + x2i;
1163
     _TEMP[ j1     ] := x0r - x2r;
1164
     _TEMP[ j1 + 1 ] := x0i - x2i;
1165

1166
     x0r := x1r - x3i; x0i := x1i + x3r;
1167
     _TEMP[ j2     ] := wn4r * ( x0r - x0i );
1168
     _TEMP[ j2 + 1 ] := wn4r * ( x0i + x0r );
1169

1170
     x0r := x1r + x3i; x0i := x1i - x3r;
1171
     _TEMP[ j3     ] := -wn4r * ( x0r + x0i );
1172
     _TEMP[ j3 + 1 ] := -wn4r * ( x0i - x0r );
1173

1174
     x0r := _TEMP[ j0 + 2 ] + _TEMP[ j2 + 2 ]; x0i := _TEMP[ j0 + 3 ] + _TEMP[ j2 + 3 ];
1175
     x1r := _TEMP[ j0 + 2 ] - _TEMP[ j2 + 2 ]; x1i := _TEMP[ j0 + 3 ] - _TEMP[ j2 + 3 ];
1176
     x2r := _TEMP[ j1 + 2 ] + _TEMP[ j3 + 2 ]; x2i := _TEMP[ j1 + 3 ] + _TEMP[ j3 + 3 ];
1177
     x3r := _TEMP[ j1 + 2 ] - _TEMP[ j3 + 2 ]; x3i := _TEMP[ j1 + 3 ] - _TEMP[ j3 + 3 ];
1178
     _TEMP[ j0 + 2 ] := x0r + x2r;
1179
     _TEMP[ j0 + 3 ] := x0i + x2i;
1180
     _TEMP[ j1 + 2 ] := x0r - x2r;
1181
     _TEMP[ j1 + 3 ] := x0i - x2i;
1182

1183
     x0r := x1r - x3i; x0i := x1i + x3r;
1184
     _TEMP[ j2 + 2 ] := wk1i * x0r - wk1r * x0i;
1185
     _TEMP[ j2 + 3 ] := wk1i * x0i + wk1r * x0r;
1186

1187
     x0r := x1r + x3i; x0i := x1i - x3r;
1188
     _TEMP[ j3 + 2 ] := wk3i * x0r + wk3r * x0i;
1189
     _TEMP[ j3 + 3 ] := wk3i * x0i - wk3r * x0r;
1190
end;
1191

1192
procedure TDiscreteTrans.cftb1st( var w:array of Double );
1193
var
1194
   j, j0, j1, j2, j3, k, m, mh :Integer;
1195
   wn4r, csc1, csc3, wk1r, wk1i, wk3r, wk3i,
1196
   wd1r, wd1i, wd3r, wd3i,
1197
   x0r, x0i, x1r, x1i, x2r, x2i, x3r, x3i,
1198
   y0r, y0i, y1r, y1i, y2r, y2i, y3r, y3i :Double;
1199
begin
1200
     mh := _TempN shr 3;
1201
     m := 2 * mh;
1202

1203
     j1 := m;
1204
     j2 := j1 + m;
1205
     j3 := j2 + m;
1206

1207
     x0r :=  _TEMP[      0 ] + _TEMP[ j2     ]; x0i := -_TEMP[      1 ] - _TEMP[ j2 + 1 ];
1208
     x1r :=  _TEMP[      0 ] - _TEMP[ j2     ]; x1i := -_TEMP[      1 ] + _TEMP[ j2 + 1 ];
1209
     x2r :=  _TEMP[ j1     ] + _TEMP[ j3     ]; x2i :=  _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1210
     x3r :=  _TEMP[ j1     ] - _TEMP[ j3     ]; x3i :=  _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1211
     _TEMP[      0 ] := x0r + x2r;
1212
     _TEMP[      1 ] := x0i - x2i;
1213
     _TEMP[ j1     ] := x0r - x2r;
1214
     _TEMP[ j1 + 1 ] := x0i + x2i;
1215
     _TEMP[ j2     ] := x1r + x3i;
1216
     _TEMP[ j2 + 1 ] := x1i + x3r;
1217
     _TEMP[ j3     ] := x1r - x3i;
1218
     _TEMP[ j3 + 1 ] := x1i - x3r;
1219

1220
     wn4r := w[ 1 ];
1221
     csc1 := w[ 2 ];
1222
     csc3 := w[ 3 ];
1223
     wd1r := 1;
1224
     wd1i := 0;
1225
     wd3r := 1;
1226
     wd3i := 0;
1227
     k := 0;
1228
     j := 2;
1229
     while j < mh - 2 do
1230
     begin
1231
          Inc( k, 4 );
1232

1233
          wk1r := csc1 * ( wd1r + w[ k     ] ); wk1i := csc1 * ( wd1i + w[ k + 1 ] );
1234
          wk3r := csc3 * ( wd3r + w[ k + 2 ] ); wk3i := csc3 * ( wd3i + w[ k + 3 ] );
1235
          wd1r := w[ k     ]; wd1i := w[ k + 1 ];
1236
          wd3r := w[ k + 2 ]; wd3i := w[ k + 3 ];
1237

1238
          j1 := j  + m;
1239
          j2 := j1 + m;
1240
          j3 := j2 + m;
1241

1242
          x0r :=  _TEMP[ j      ] + _TEMP[ j2     ];
1243
          x0i := -_TEMP[ j  + 1 ] - _TEMP[ j2 + 1 ];
1244
          x1r :=  _TEMP[ j      ] - _TEMP[ j2     ];
1245
          x1i := -_TEMP[ j  + 1 ] + _TEMP[ j2 + 1 ];
1246
          y0r :=  _TEMP[ j  + 2 ] + _TEMP[ j2 + 2 ];
1247
          y0i := -_TEMP[ j  + 3 ] - _TEMP[ j2 + 3 ];
1248
          y1r :=  _TEMP[ j  + 2 ] - _TEMP[ j2 + 2 ];
1249
          y1i := -_TEMP[ j  + 3 ] + _TEMP[ j2 + 3 ];
1250
          x2r :=  _TEMP[ j1     ] + _TEMP[ j3     ];
1251
          x2i :=  _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1252
          x3r :=  _TEMP[ j1     ] - _TEMP[ j3     ];
1253
          x3i :=  _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1254
          y2r :=  _TEMP[ j1 + 2 ] + _TEMP[ j3 + 2 ];
1255
          y2i :=  _TEMP[ j1 + 3 ] + _TEMP[ j3 + 3 ];
1256
          y3r :=  _TEMP[ j1 + 2 ] - _TEMP[ j3 + 2 ];
1257
          y3i :=  _TEMP[ j1 + 3 ] - _TEMP[ j3 + 3 ];
1258
          _TEMP[ j      ] := x0r + x2r;
1259
          _TEMP[ j  + 1 ] := x0i - x2i;
1260
          _TEMP[ j  + 2 ] := y0r + y2r;
1261
          _TEMP[ j  + 3 ] := y0i - y2i;
1262
          _TEMP[ j1     ] := x0r - x2r;
1263
          _TEMP[ j1 + 1 ] := x0i + x2i;
1264
          _TEMP[ j1 + 2 ] := y0r - y2r;
1265
          _TEMP[ j1 + 3 ] := y0i + y2i;
1266

1267
          x0r := x1r + x3i; x0i := x1i + x3r;
1268
          _TEMP[ j2     ] := wk1r * x0r - wk1i * x0i;
1269
          _TEMP[ j2 + 1 ] := wk1r * x0i + wk1i * x0r;
1270

1271
          x0r := y1r + y3i; x0i := y1i + y3r;
1272
          _TEMP[ j2 + 2 ] := wd1r * x0r - wd1i * x0i;
1273
          _TEMP[ j2 + 3 ] := wd1r * x0i + wd1i * x0r;
1274

1275
          x0r := x1r - x3i; x0i := x1i - x3r;
1276
          _TEMP[ j3     ] := wk3r * x0r + wk3i * x0i;
1277
          _TEMP[ j3 + 1 ] := wk3r * x0i - wk3i * x0r;
1278

1279
          x0r := y1r - y3i; x0i := y1i - y3r;
1280
          _TEMP[ j3 + 2 ] := wd3r * x0r + wd3i * x0i;
1281
          _TEMP[ j3 + 3 ] := wd3r * x0i - wd3i * x0r;
1282

1283
          j0 :=  m - j;
1284
          j1 := j0 + m;
1285
          j2 := j1 + m;
1286
          j3 := j2 + m;
1287

1288
          x0r :=  _TEMP[ j0     ] + _TEMP[ j2     ]; x0i := -_TEMP[ j0 + 1 ] - _TEMP[ j2 + 1 ];
1289
          x1r :=  _TEMP[ j0     ] - _TEMP[ j2     ]; x1i := -_TEMP[ j0 + 1 ] + _TEMP[ j2 + 1 ];
1290
          y0r :=  _TEMP[ j0 - 2 ] + _TEMP[ j2 - 2 ]; y0i := -_TEMP[ j0 - 1 ] - _TEMP[ j2 - 1 ];
1291
          y1r :=  _TEMP[ j0 - 2 ] - _TEMP[ j2 - 2 ]; y1i := -_TEMP[ j0 - 1 ] + _TEMP[ j2 - 1 ];
1292
          x2r :=  _TEMP[ j1     ] + _TEMP[ j3     ]; x2i :=  _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1293
          x3r :=  _TEMP[ j1     ] - _TEMP[ j3     ]; x3i :=  _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1294
          y2r :=  _TEMP[ j1 - 2 ] + _TEMP[ j3 - 2 ]; y2i :=  _TEMP[ j1 - 1 ] + _TEMP[ j3 - 1 ];
1295
          y3r :=  _TEMP[ j1 - 2 ] - _TEMP[ j3 - 2 ]; y3i :=  _TEMP[ j1 - 1 ] - _TEMP[ j3 - 1 ];
1296
          _TEMP[ j0     ] := x0r + x2r;
1297
          _TEMP[ j0 + 1 ] := x0i - x2i;
1298
          _TEMP[ j0 - 2 ] := y0r + y2r;
1299
          _TEMP[ j0 - 1 ] := y0i - y2i;
1300
          _TEMP[ j1     ] := x0r - x2r;
1301
          _TEMP[ j1 + 1 ] := x0i + x2i;
1302
          _TEMP[ j1 - 2 ] := y0r - y2r;
1303
          _TEMP[ j1 - 1 ] := y0i + y2i;
1304

1305
          x0r := x1r + x3i; x0i := x1i + x3r;
1306
          _TEMP[ j2 ] := wk1i * x0r - wk1r * x0i;
1307
          _TEMP[ j2 + 1 ] := wk1i * x0i + wk1r * x0r;
1308

1309
          x0r := y1r + y3i; x0i := y1i + y3r;
1310
          _TEMP[ j2 - 2 ] := wd1i * x0r - wd1r * x0i;
1311
          _TEMP[ j2 - 1 ] := wd1i * x0i + wd1r * x0r;
1312

1313
          x0r := x1r - x3i; x0i := x1i - x3r;
1314
          _TEMP[ j3 ] := wk3i * x0r + wk3r * x0i;
1315
          _TEMP[ j3 + 1 ] := wk3i * x0i - wk3r * x0r;
1316

1317
          x0r := y1r - y3i; x0i := y1i - y3r;
1318
          _TEMP[ j3 - 2 ] := wd3i * x0r + wd3r * x0i;
1319
          _TEMP[ j3 - 1 ] := wd3i * x0i - wd3r * x0r;
1320

1321
          Inc( j, 4 )
1322
     end;
1323

1324
     wk1r := csc1 * ( wd1r + wn4r );
1325
     wk1i := csc1 * ( wd1i + wn4r );
1326
     wk3r := csc3 * ( wd3r - wn4r );
1327
     wk3i := csc3 * ( wd3i - wn4r );
1328

1329
     j0 := mh;
1330
     j1 := j0 + m;
1331
     j2 := j1 + m;
1332
     j3 := j2 + m;
1333

1334
     x0r :=  _TEMP[ j0 - 2 ] + _TEMP[ j2 - 2 ];  x0i := -_TEMP[ j0 - 1 ] - _TEMP[ j2 - 1 ];
1335
     x1r :=  _TEMP[ j0 - 2 ] - _TEMP[ j2 - 2 ];  x1i := -_TEMP[ j0 - 1 ] + _TEMP[ j2 - 1 ];
1336
     x2r :=  _TEMP[ j1 - 2 ] + _TEMP[ j3 - 2 ];  x2i :=  _TEMP[ j1 - 1 ] + _TEMP[ j3 - 1 ];
1337
     x3r :=  _TEMP[ j1 - 2 ] - _TEMP[ j3 - 2 ];  x3i :=  _TEMP[ j1 - 1 ] - _TEMP[ j3 - 1 ];
1338
     _TEMP[ j0 - 2 ] := x0r + x2r;
1339
     _TEMP[ j0 - 1 ] := x0i - x2i;
1340
     _TEMP[ j1 - 2 ] := x0r - x2r;
1341
     _TEMP[ j1 - 1 ] := x0i + x2i;
1342

1343
     x0r := x1r + x3i; x0i := x1i + x3r;
1344
     _TEMP[ j2 - 2 ] := wk1r * x0r - wk1i * x0i;
1345
     _TEMP[ j2 - 1 ] := wk1r * x0i + wk1i * x0r;
1346

1347
     x0r := x1r - x3i; x0i := x1i - x3r;
1348
     _TEMP[ j3 - 2 ] := wk3r * x0r + wk3i * x0i;
1349
     _TEMP[ j3 - 1 ] := wk3r * x0i - wk3i * x0r;
1350

1351
     x0r :=  _TEMP[ j0     ] + _TEMP[ j2     ];  x0i := -_TEMP[ j0 + 1 ] - _TEMP[ j2 + 1 ];
1352
     x1r :=  _TEMP[ j0     ] - _TEMP[ j2     ];  x1i := -_TEMP[ j0 + 1 ] + _TEMP[ j2 + 1 ];
1353
     x2r :=  _TEMP[ j1     ] + _TEMP[ j3     ];  x2i :=  _TEMP[ j1 + 1 ] + _TEMP[ j3 + 1 ];
1354
     x3r :=  _TEMP[ j1     ] - _TEMP[ j3     ];  x3i :=  _TEMP[ j1 + 1 ] - _TEMP[ j3 + 1 ];
1355
     _TEMP[ j0     ] := x0r + x2r;
1356
     _TEMP[ j0 + 1 ] := x0i - x2i;
1357
     _TEMP[ j1     ] := x0r - x2r;
1358
     _TEMP[ j1 + 1 ] := x0i + x2i;
1359

1360
     x0r := x1r + x3i; x0i := x1i + x3r;
1361
     _TEMP[ j2 ] := wn4r * ( x0r - x0i );
1362
     _TEMP[ j2 + 1 ] := wn4r * ( x0i + x0r );
1363

1364
     x0r := x1r - x3i; x0i := x1i - x3r;
1365
     _TEMP[ j3 ] := -wn4r * ( x0r + x0i );
1366
     _TEMP[ j3 + 1 ] := -wn4r * ( x0i - x0r );
1367

1368
     x0r :=  _TEMP[ j0 + 2 ] + _TEMP[ j2 + 2 ];  x0i := -_TEMP[ j0 + 3 ] - _TEMP[ j2 + 3 ];
1369
     x1r :=  _TEMP[ j0 + 2 ] - _TEMP[ j2 + 2 ];  x1i := -_TEMP[ j0 + 3 ] + _TEMP[ j2 + 3 ];
1370
     x2r :=  _TEMP[ j1 + 2 ] + _TEMP[ j3 + 2 ];  x2i :=  _TEMP[ j1 + 3 ] + _TEMP[ j3 + 3 ];
1371
     x3r :=  _TEMP[ j1 + 2 ] - _TEMP[ j3 + 2 ];  x3i :=  _TEMP[ j1 + 3 ] - _TEMP[ j3 + 3 ];
1372
     _TEMP[ j0 + 2 ] := x0r + x2r;
1373
     _TEMP[ j0 + 3 ] := x0i - x2i;
1374
     _TEMP[ j1 + 2 ] := x0r - x2r;
1375
     _TEMP[ j1 + 3 ] := x0i + x2i;
1376

1377
     x0r := x1r + x3i; x0i := x1i + x3r;
1378
     _TEMP[ j2 + 2 ] := wk1i * x0r - wk1r * x0i;
1379
     _TEMP[ j2 + 3 ] := wk1i * x0i + wk1r * x0r;
1380

1381
     x0r := x1r - x3i; x0i := x1i - x3r;
1382
     _TEMP[ j3 + 2 ] := wk3i * x0r + wk3r * x0i;
1383
     _TEMP[ j3 + 3 ] := wk3i * x0i - wk3r * x0r;
1384
end;
1385

1386
////////////////////////////////////////////////////////////////////////////////////////////////////
1387

1388
procedure TDiscreteTrans.cftmdl1( const n:Integer; var a:array of Double; var w:array of Double );
1389
var
1390
   j, j0, j1, j2, j3, k, m, mh :Integer;
1391
   wn4r, wk1r, wk1i, wk3r, wk3i,
1392
   x0r, x0i, x1r, x1i, x2r, x2i, x3r, x3i :Double;
1393
begin
1394
    mh := n shr 3;
1395
    m := 2 * mh;
1396

1397
    j1 :=      m;
1398
    j2 := j1 + m;
1399
    j3 := j2 + m;
1400

1401
    x0r := a[      0 ] + a[ j2     ];
1402
    x0i := a[      1 ] + a[ j2 + 1 ];
1403
    x1r := a[      0 ] - a[ j2     ];
1404
    x1i := a[      1 ] - a[ j2 + 1 ];
1405
    x2r := a[ j1     ] + a[ j3     ];
1406
    x2i := a[ j1 + 1 ] + a[ j3 + 1 ];
1407
    x3r := a[ j1     ] - a[ j3     ];
1408
    x3i := a[ j1 + 1 ] - a[ j3 + 1 ];
1409
    a[      0 ] := x0r + x2r;
1410
    a[      1 ] := x0i + x2i;
1411
    a[ j1     ] := x0r - x2r;
1412
    a[ j1 + 1 ] := x0i - x2i;
1413
    a[ j2     ] := x1r - x3i;
1414
    a[ j2 + 1 ] := x1i + x3r;
1415
    a[ j3     ] := x1r + x3i;
1416
    a[ j3 + 1 ] := x1i - x3r;
1417

1418
    wn4r := w[ 1 ];
1419
    k := 0;
1420
    j := 2;
1421
    while j < mh do
1422
    begin
1423
        Inc( k, 4 );
1424

1425
        wk1r := w[ k     ];
1426
        wk1i := w[ k + 1 ];
1427
        wk3r := w[ k + 2 ];
1428
        wk3i := w[ k + 3 ];
1429

1430
        j1 := j  + m;
1431
        j2 := j1 + m;
1432
        j3 := j2 + m;
1433

1434
        x0r := a[ j      ] + a[ j2     ];
1435
        x0i := a[ j  + 1 ] + a[ j2 + 1 ];
1436
        x1r := a[ j      ] - a[ j2     ];
1437
        x1i := a[ j  + 1 ] - a[ j2 + 1 ];
1438
        x2r := a[ j1     ] + a[ j3     ];
1439
        x2i := a[ j1 + 1 ] + a[ j3 + 1 ];
1440
        x3r := a[ j1     ] - a[ j3     ];
1441
        x3i := a[ j1 + 1 ] - a[ j3 + 1 ];
1442
        a[ j      ] := x0r + x2r;
1443
        a[ j  + 1 ] := x0i + x2i;
1444
        a[ j1     ] := x0r - x2r;
1445
        a[ j1 + 1 ] := x0i - x2i;
1446

1447
        x0r := x1r - x3i;
1448
        x0i := x1i + x3r;
1449
        a[ j2     ] := wk1r * x0r - wk1i * x0i;
1450
        a[ j2 + 1 ] := wk1r * x0i + wk1i * x0r;
1451

1452
        x0r := x1r + x3i;
1453
        x0i := x1i - x3r;
1454
        a[ j3     ] := wk3r * x0r + wk3i * x0i;
1455
        a[ j3 + 1 ] := wk3r * x0i - wk3i * x0r;
1456

1457
        j0 := m  - j;
1458
        j1 := j0 + m;
1459
        j2 := j1 + m;
1460
        j3 := j2 + m;
1461

1462
        x0r := a[ j0     ] + a[ j2     ];
1463
        x0i := a[ j0 + 1 ] + a[ j2 + 1 ];
1464
        x1r := a[ j0     ] - a[ j2     ];
1465
        x1i := a[ j0 + 1 ] - a[ j2 + 1 ];
1466
        x2r := a[ j1     ] + a[ j3     ];
1467
        x2i := a[ j1 + 1 ] + a[ j3 + 1 ];
1468
        x3r := a[ j1     ] - a[ j3     ];
1469
        x3i := a[ j1 + 1 ] - a[ j3 + 1 ];
1470
        a[ j0     ] := x0r + x2r;
1471
        a[ j0 + 1 ] := x0i + x2i;
1472
        a[ j1     ] := x0r - x2r;
1473
        a[ j1 + 1 ] := x0i - x2i;
1474

1475
        x0r := x1r - x3i;
1476
        x0i := x1i + x3r;
1477
        a[ j2     ] := wk1i * x0r - wk1r * x0i;
1478
        a[ j2 + 1 ] := wk1i * x0i + wk1r * x0r;
1479

1480
        x0r := x1r + x3i;
1481
        x0i := x1i - x3r;
1482
        a[ j3     ] := wk3i * x0r + wk3r * x0i;
1483
        a[ j3 + 1 ] := wk3i * x0i - wk3r * x0r;
1484

1485
        Inc( j, 2 )
1486
    end;
1487

1488
    j0 := mh;
1489
    j1 := j0 + m;
1490
    j2 := j1 + m;
1491
    j3 := j2 + m;
1492

1493
    x0r := a[ j0     ] + a[ j2     ];
1494
    x0i := a[ j0 + 1 ] + a[ j2 + 1 ];
1495
    x1r := a[ j0     ] - a[ j2     ];
1496
    x1i := a[ j0 + 1 ] - a[ j2 + 1 ];
1497
    x2r := a[ j1     ] + a[ j3     ];
1498
    x2i := a[ j1 + 1 ] + a[ j3 + 1 ];
1499
    x3r := a[ j1     ] - a[ j3     ];
1500
    x3i := a[ j1 + 1 ] - a[ j3 + 1 ];
1501
    a[ j0     ] := x0r + x2r;
1502
    a[ j0 + 1 ] := x0i + x2i;
1503
    a[ j1     ] := x0r - x2r;
1504
    a[ j1 + 1 ] := x0i - x2i;
1505

1506
    x0r := x1r - x3i;
1507
    x0i := x1i + x3r;
1508
    a[ j2     ] := +wn4r * ( x0r - x0i );
1509
    a[ j2 + 1 ] := +wn4r * ( x0i + x0r );
1510

1511
    x0r := x1r + x3i;
1512
    x0i := x1i - x3r;
1513
    a[ j3     ] := -wn4r * ( x0r + x0i );
1514
    a[ j3 + 1 ] := -wn4r * ( x0i - x0r );
1515
end;
1516

1517
procedure TDiscreteTrans.cftmdl2( const n:Integer; var a:array of Double; var w:array of Double );
1518
var
1519
   j, j0, j1, j2, j3, k, kr, m, mh :Integer;
1520
   wn4r, wk1r, wk1i, wk3r, wk3i, wd1r, wd1i, wd3r, wd3i,
1521
   x0r, x0i, x1r, x1i, x2r, x2i, x3r, x3i, y0r, y0i, y2r, y2i :Double;
1522
 begin
1523
    mh := n shr 3;
1524
    m := 2 * mh;
1525

1526
    wn4r := w[ 1 ];
1527

1528
    j1 := m;
1529
    j2 := j1 + m;
1530
    j3 := j2 + m;
1531

1532
    x0r := a[      0 ] - a[ j2 + 1 ];
1533
    x0i := a[      1 ] + a[ j2     ];
1534
    x1r := a[      0 ] + a[ j2 + 1 ];
1535
    x1i := a[      1 ] - a[ j2     ];
1536
    x2r := a[ j1     ] - a[ j3 + 1 ];
1537
    x2i := a[ j1 + 1 ] + a[ j3     ];
1538
    x3r := a[ j1     ] + a[ j3 + 1 ];
1539
    x3i := a[ j1 + 1 ] - a[ j3     ];
1540
    y0r := wn4r * ( x2r - x2i );
1541
    y0i := wn4r * ( x2i + x2r );
1542
    a[      0 ] := x0r + y0r;
1543
    a[      1 ] := x0i + y0i;
1544
    a[ j1     ] := x0r - y0r;
1545
    a[ j1 + 1 ] := x0i - y0i;
1546

1547
    y0r := wn4r * ( x3r - x3i );
1548
    y0i := wn4r * ( x3i + x3r );
1549
    a[ j2     ] := x1r - y0i;
1550
    a[ j2 + 1 ] := x1i + y0r;
1551
    a[ j3     ] := x1r + y0i;
1552
    a[ j3 + 1 ] := x1i - y0r;
1553

1554
    k := 0;
1555
    kr := 2 * m;
1556
    j := 2;
1557
    while j < mh do
1558
    begin
1559
        Inc( k, 4 );
1560

1561
        wk1r := w[ k     ];
1562
        wk1i := w[ k + 1 ];
1563
        wk3r := w[ k + 2 ];
1564
        wk3i := w[ k + 3 ];
1565

1566
        Dec( kr, 4 );
1567
        wd1i := w[ kr     ];
1568
        wd1r := w[ kr + 1 ];
1569
        wd3i := w[ kr + 2 ];
1570
        wd3r := w[ kr + 3 ];
1571

1572
        j1 := j  + m;
1573
        j2 := j1 + m;
1574
        j3 := j2 + m;
1575

1576
        x0r := a[ j      ] - a[ j2 + 1 ];
1577
        x0i := a[ j  + 1 ] + a[ j2     ];
1578
        x1r := a[ j      ] + a[ j2 + 1 ];
1579
        x1i := a[ j  + 1 ] - a[ j2     ];
1580
        x2r := a[ j1     ] - a[ j3 + 1 ];
1581
        x2i := a[ j1 + 1 ] + a[ j3     ];
1582
        x3r := a[ j1     ] + a[ j3 + 1 ];
1583
        x3i := a[ j1 + 1 ] - a[ j3     ];
1584

1585
        y0r := wk1r * x0r - wk1i * x0i;
1586
        y0i := wk1r * x0i + wk1i * x0r;
1587
        y2r := wd1r * x2r - wd1i * x2i;
1588
        y2i := wd1r * x2i + wd1i * x2r;
1589
        a[ j      ] := y0r + y2r;
1590
        a[ j  + 1 ] := y0i + y2i;
1591
        a[ j1     ] := y0r - y2r;
1592
        a[ j1 + 1 ] := y0i - y2i;
1593

1594
        y0r := wk3r * x1r + wk3i * x1i;
1595
        y0i := wk3r * x1i - wk3i * x1r;
1596
        y2r := wd3r * x3r + wd3i * x3i;
1597
        y2i := wd3r * x3i - wd3i * x3r;
1598
        a[ j2     ] := y0r + y2r;
1599
        a[ j2 + 1 ] := y0i + y2i;
1600
        a[ j3     ] := y0r - y2r;
1601
        a[ j3 + 1 ] := y0i - y2i;
1602

1603
        j0 := m  - j;
1604
        j1 := j0 + m;
1605
        j2 := j1 + m;
1606
        j3 := j2 + m;
1607

1608
        x0r := a[ j0     ] - a[ j2 + 1 ];
1609
        x0i := a[ j0 + 1 ] + a[ j2     ];
1610
        x1r := a[ j0     ] + a[ j2 + 1 ];
1611
        x1i := a[ j0 + 1 ] - a[ j2     ];
1612
        x2r := a[ j1     ] - a[ j3 + 1 ];
1613
        x2i := a[ j1 + 1 ] + a[ j3     ];
1614
        x3r := a[ j1     ] + a[ j3 + 1 ];
1615
        x3i := a[ j1 + 1 ] - a[ j3     ];
1616

1617
        y0r := wd1i * x0r - wd1r * x0i;
1618
        y0i := wd1i * x0i + wd1r * x0r;
1619
        y2r := wk1i * x2r - wk1r * x2i;
1620
        y2i := wk1i * x2i + wk1r * x2r;
1621
        a[ j0     ] := y0r + y2r;
1622
        a[ j0 + 1 ] := y0i + y2i;
1623
        a[ j1     ] := y0r - y2r;
1624
        a[ j1 + 1 ] := y0i - y2i;
1625

1626
        y0r := wd3i * x1r + wd3r * x1i;
1627
        y0i := wd3i * x1i - wd3r * x1r;
1628
        y2r := wk3i * x3r + wk3r * x3i;
1629
        y2i := wk3i * x3i - wk3r * x3r;
1630
        a[ j2     ] := y0r + y2r;
1631
        a[ j2 + 1 ] := y0i + y2i;
1632
        a[ j3     ] := y0r - y2r;
1633
        a[ j3 + 1 ] := y0i - y2i;
1634

1635
        Inc( j, 2 )
1636
    end;
1637

1638
    wk1r := w[ m ];
1639
    wk1i := w[ m + 1 ];
1640

1641
    j0 := mh;
1642
    j1 := j0 + m;
1643
    j2 := j1 + m;
1644
    j3 := j2 + m;
1645

1646
    x0r := a[ j0     ] - a[ j2 + 1 ];
1647
    x0i := a[ j0 + 1 ] + a[ j2     ];
1648
    x1r := a[ j0     ] + a[ j2 + 1 ];
1649
    x1i := a[ j0 + 1 ] - a[ j2     ];
1650
    x2r := a[ j1     ] - a[ j3 + 1 ];
1651
    x2i := a[ j1 + 1 ] + a[ j3     ];
1652
    x3r := a[ j1     ] + a[ j3 + 1 ];
1653
    x3i := a[ j1 + 1 ] - a[ j3     ];
1654
    y0r := wk1r * x0r - wk1i * x0i;
1655
    y0i := wk1r * x0i + wk1i * x0r;
1656
    y2r := wk1i * x2r - wk1r * x2i;
1657
    y2i := wk1i * x2i + wk1r * x2r;
1658
    a[ j0     ] := y0r + y2r;
1659
    a[ j0 + 1 ] := y0i + y2i;
1660
    a[ j1     ] := y0r - y2r;
1661
    a[ j1 + 1 ] := y0i - y2i;
1662

1663
    y0r := wk1i * x1r - wk1r * x1i; y0i := wk1i * x1i + wk1r * x1r;
1664
    y2r := wk1r * x3r - wk1i * x3i; y2i := wk1r * x3i + wk1i * x3r;
1665
    a[ j2     ] := y0r - y2r;
1666
    a[ j2 + 1 ] := y0i - y2i;
1667
    a[ j3     ] := y0r + y2r;
1668
    a[ j3 + 1 ] := y0i + y2i;
1669
end;
1670

1671
////////////////////////////////////////////////////////////////////////////////////////////////////
1672

1673
function TDiscreteTrans.cfttree( const n,j,k:Integer ) :Integer;
1674
var
1675
   i, isplt, m :Integer;
1676
begin
1677
     if k and 3 <> 0 then
1678
     begin
1679
          isplt := k and 1;
1680

1681
          if isplt <> 0 then cftmdl1( n, _TEMP[ j - n ], _W[ _NW - ( n shr 1 ) ] )
1682
                        else cftmdl2( n, _TEMP[ j - n ], _W[ _NW -   n         ] )
1683
     end
1684
     else
1685
     begin
1686
          m := n;
1687
          i := k;
1688
          while i and 3 = 0 do
1689
          begin
1690
               m := m shl 2;
1691
               i := i shr 2;
1692
          end;
1693

1694
          isplt := i and 1;
1695

1696
          if isplt <> 0 then
1697
          begin
1698
               while m > 128 do
1699
               begin
1700
                    cftmdl1( m, _TEMP[ j - m ], _W[ _NW - ( m shr 1 ) ]);
1701

1702
                    m := m shr 2
1703
               end
1704
          end
1705
          else
1706
          begin
1707
               while m > 128 do
1708
               begin
1709
                    cftmdl2( m, _TEMP[ j - m ], _W[ _NW - m ] );
1710

1711
                    m := m shr 2
1712
               end
1713
          end
1714
     end;
1715

1716
     Result := isplt
1717
end;
1718

1719

1720
////////////////////////////////////////////////////////////////////////////////////////////////////
1721

1722
procedure TDiscreteTrans.cftf161( var a:array of Double; var w:array of Double );
1723
var
1724
   wn4r, wk1r, wk1i,
1725
   x00r, x00i, x01r, x01i, x02r, x02i, x03r, x03i,
1726
   y00r, y00i, y01r, y01i, y02r, y02i, y03r, y03i,
1727
   y04r, y04i, y05r, y05i, y06r, y06i, y07r, y07i,
1728
   y08r, y08i, y09r, y09i, y10r, y10i, y11r, y11i,
1729
   y12r, y12i, y13r, y13i, y14r, y14i, y15r, y15i :Double;
1730
begin
1731
     wn4r := w[ 1 ];
1732
     wk1r := w[ 2 ];
1733
     wk1i := w[ 3 ];
1734

1735
     x00r := a[ 00 ] + a[ 16 ]; x00i := a[ 01 ] + a[ 17 ];
1736
     x01r := a[ 00 ] - a[ 16 ]; x01i := a[ 01 ] - a[ 17 ];
1737
     x02r := a[ 08 ] + a[ 24 ]; x02i := a[ 09 ] + a[ 25 ];
1738
     x03r := a[ 08 ] - a[ 24 ]; x03i := a[ 09 ] - a[ 25 ];
1739
     y00r := x00r + x02r;       y00i := x00i + x02i;
1740
     y04r := x00r - x02r;       y04i := x00i - x02i;
1741
     y08r := x01r - x03i;       y08i := x01i + x03r;
1742
     y12r := x01r + x03i;       y12i := x01i - x03r;
1743
     x00r := a[ 02 ] + a[ 18 ];
1744
     x00i := a[ 03 ] + a[ 19 ];
1745
     x01r := a[ 02 ] - a[ 18 ];
1746
     x01i := a[ 03 ] - a[ 19 ];
1747
     x02r := a[ 10 ] + a[ 26 ];
1748
     x02i := a[ 11 ] + a[ 27 ];
1749
     x03r := a[ 10 ] - a[ 26 ];
1750
     x03i := a[ 11 ] - a[ 27 ];
1751

1752
     y01r := x00r + x02r;               y01i := x00i + x02i;
1753
     y05r := x00r - x02r;               y05i := x00i - x02i;
1754
     x00r := x01r - x03i;               x00i := x01i + x03r;
1755
     y09r := wk1r * x00r - wk1i * x00i; y09i := wk1r * x00i + wk1i * x00r;
1756
     x00r := x01r + x03i;               x00i := x01i - x03r;
1757
     y13r := wk1i * x00r - wk1r * x00i; y13i := wk1i * x00i + wk1r * x00r;
1758
     x00r := a[ 4 ] + a[ 20 ];
1759
     x00i := a[ 5 ] + a[ 21 ];
1760
     x01r := a[ 4 ] - a[ 20 ];
1761
     x01i := a[ 5 ] - a[ 21 ];
1762
     x02r := a[ 12 ] + a[ 28 ];
1763
     x02i := a[ 13 ] + a[ 29 ];
1764
     x03r := a[ 12 ] - a[ 28 ];
1765
     x03i := a[ 13 ] - a[ 29 ];
1766

1767
     y02r := x00r + x02r;              y02i := x00i + x02i;
1768
     y06r := x00r - x02r;              y06i := x00i - x02i;
1769
     x00r := x01r - x03i;              x00i := x01i + x03r;
1770
     y10r := wn4r * (  x00r - x00i  ); y10i := wn4r * (  x00i + x00r  );
1771
     x00r := x01r + x03i;              x00i := x01i - x03r;
1772
     y14r := wn4r * (  x00r + x00i  ); y14i := wn4r * (  x00i - x00r  );
1773
     x00r := a[ 06 ] + a[ 22 ];
1774
     x00i := a[ 07 ] + a[ 23 ];
1775
     x01r := a[ 06 ] - a[ 22 ];
1776
     x01i := a[ 07 ] - a[ 23 ];
1777
     x02r := a[ 14 ] + a[ 30 ];
1778
     x02i := a[ 15 ] + a[ 31 ];
1779
     x03r := a[ 14 ] - a[ 30 ];
1780
     x03i := a[ 15 ] - a[ 31 ];
1781

1782
     y03r := x00r + x02r;               y03i := x00i + x02i;
1783
     y07r := x00r - x02r;               y07i := x00i - x02i;
1784
     x00r := x01r - x03i;               x00i := x01i + x03r;
1785
     y11r := wk1i * x00r - wk1r * x00i; y11i := wk1i * x00i + wk1r * x00r;
1786
     x00r := x01r + x03i;               x00i := x01i - x03r;
1787
     y15r := wk1r * x00r - wk1i * x00i; y15i := wk1r * x00i + wk1i * x00r;
1788
     x00r := y12r - y14r;               x00i := y12i - y14i;
1789
     x01r := y12r + y14r;               x01i := y12i + y14i;
1790
     x02r := y13r - y15r;               x02i := y13i - y15i;
1791
     x03r := y13r + y15r;               x03i := y13i + y15i;
1792
     a[ 24 ] := x00r + x02r;
1793
     a[ 25 ] := x00i + x02i;
1794
     a[ 26 ] := x00r - x02r;
1795
     a[ 27 ] := x00i - x02i;
1796
     a[ 28 ] := x01r - x03i;
1797
     a[ 29 ] := x01i + x03r;
1798
     a[ 30 ] := x01r + x03i;
1799
     a[ 31 ] := x01i - x03r;
1800

1801
     x00r := y08r + y10r; x00i := y08i + y10i;
1802
     x01r := y08r - y10r; x01i := y08i - y10i;
1803
     x02r := y09r + y11r; x02i := y09i + y11i;
1804
     x03r := y09r - y11r; x03i := y09i - y11i;
1805
     a[ 16 ] := x00r + x02r;
1806
     a[ 17 ] := x00i + x02i;
1807
     a[ 18 ] := x00r - x02r;
1808
     a[ 19 ] := x00i - x02i;
1809
     a[ 20 ] := x01r - x03i;
1810
     a[ 21 ] := x01i + x03r;
1811
     a[ 22 ] := x01r + x03i;
1812
     a[ 23 ] := x01i - x03r;
1813

1814
     x00r := y05r - y07i;            x00i := y05i + y07r;
1815
     x02r := wn4r * ( x00r - x00i ); x02i := wn4r * ( x00i + x00r );
1816
     x00r := y05r + y07i;            x00i := y05i - y07r;
1817
     x03r := wn4r * ( x00r - x00i ); x03i := wn4r * ( x00i + x00r );
1818
     x00r := y04r - y06i;            x00i := y04i + y06r;
1819
     x01r := y04r + y06i;            x01i := y04i - y06r;
1820
     a[ 08 ] := x00r + x02r;
1821
     a[ 09 ] := x00i + x02i;
1822
     a[ 10 ] := x00r - x02r;
1823
     a[ 11 ] := x00i - x02i;
1824
     a[ 12 ] := x01r - x03i;
1825
     a[ 13 ] := x01i + x03r;
1826
     a[ 14 ] := x01r + x03i;
1827
     a[ 15 ] := x01i - x03r;
1828

1829
     x00r := y00r + y02r; x00i := y00i + y02i;
1830
     x01r := y00r - y02r; x01i := y00i - y02i;
1831
     x02r := y01r + y03r; x02i := y01i + y03i;
1832
     x03r := y01r - y03r; x03i := y01i - y03i;
1833
     a[ 00 ] := x00r + x02r;
1834
     a[ 01 ] := x00i + x02i;
1835
     a[ 02 ] := x00r - x02r;
1836
     a[ 03 ] := x00i - x02i;
1837
     a[ 04 ] := x01r - x03i;
1838
     a[ 05 ] := x01i + x03r;
1839
     a[ 06 ] := x01r + x03i;
1840
     a[ 07 ] := x01i - x03r;
1841
end;
1842

1843
procedure TDiscreteTrans.cftf162( var a:array of Double; var w:array of Double );
1844
var
1845
   wn4r, wk1r, wk1i, wk2r, wk2i, wk3r, wk3i,
1846
   x0r, x0i, x1r, x1i, x2r, x2i,
1847
   y0r, y0i, y1r, y1i, y2r, y2i, y3r, y3i,
1848
   y4r, y4i, y5r, y5i, y6r, y6i, y7r, y7i,
1849
   y8r, y8i, y9r, y9i, y10r, y10i, y11r, y11i,
1850
   y12r, y12i, y13r, y13i, y14r, y14i, y15r, y15i :Double;
1851
begin
1852
     wn4r := +w[ 1 ];
1853
     wk1r := +w[ 4 ];
1854
     wk1i := +w[ 5 ];
1855
     wk3r := +w[ 6 ];
1856
     wk3i := -w[ 7 ];
1857
     wk2r := +w[ 8 ];
1858
     wk2i := +w[ 9 ];
1859

1860
    x1r := a[ 0 ] - a[ 17 ];
1861
    x1i := a[ 1 ] + a[ 16 ];
1862
    x0r := a[ 8 ] - a[ 25 ];
1863
    x0i := a[ 9 ] + a[ 24 ];
1864
    x2r := wn4r * ( x0r - x0i );
1865
    x2i := wn4r * ( x0i + x0r );
1866
    y0r := x1r + x2r;
1867
    y0i := x1i + x2i;
1868
    y4r := x1r - x2r;
1869
    y4i := x1i - x2i;
1870
    x1r := a[ 0 ] + a[ 17 ];
1871
    x1i := a[ 1 ] - a[ 16 ];
1872
    x0r := a[ 8 ] + a[ 25 ];
1873
    x0i := a[ 9 ] - a[ 24 ];
1874
    x2r := wn4r * ( x0r - x0i );
1875
    x2i := wn4r * ( x0i + x0r );
1876
    y8r := x1r - x2i;
1877
    y8i := x1i + x2r;
1878
    y12r := x1r + x2i;
1879
    y12i := x1i - x2r;
1880
    x0r := a[ 2 ] - a[ 19 ];
1881
    x0i := a[ 3 ] + a[ 18 ];
1882
    x1r := wk1r * x0r - wk1i * x0i;
1883
    x1i := wk1r * x0i + wk1i * x0r;
1884
    x0r := a[ 10 ] - a[ 27 ];
1885
    x0i := a[ 11 ] + a[ 26 ];
1886
    x2r := wk3i * x0r - wk3r * x0i;
1887
    x2i := wk3i * x0i + wk3r * x0r;
1888
    y1r := x1r + x2r;
1889
    y1i := x1i + x2i;
1890
    y5r := x1r - x2r;
1891
    y5i := x1i - x2i;
1892
    x0r := a[ 2 ] + a[ 19 ];
1893
    x0i := a[ 3 ] - a[ 18 ];
1894
    x1r := wk3r * x0r - wk3i * x0i;
1895
    x1i := wk3r * x0i + wk3i * x0r;
1896
    x0r := a[ 10 ] + a[ 27 ];
1897
    x0i := a[ 11 ] - a[ 26 ];
1898
    x2r := wk1r * x0r + wk1i * x0i;
1899
    x2i := wk1r * x0i - wk1i * x0r;
1900
    y9r := x1r - x2r;
1901
    y9i := x1i - x2i;
1902
    y13r := x1r + x2r;
1903
    y13i := x1i + x2i;
1904
    x0r := a[ 4 ] - a[ 21 ];
1905
    x0i := a[ 5 ] + a[ 20 ];
1906
    x1r := wk2r * x0r - wk2i * x0i;
1907
    x1i := wk2r * x0i + wk2i * x0r;
1908
    x0r := a[ 12 ] - a[ 29 ];
1909
    x0i := a[ 13 ] + a[ 28 ];
1910
    x2r := wk2i * x0r - wk2r * x0i;
1911
    x2i := wk2i * x0i + wk2r * x0r;
1912
    y2r := x1r + x2r;
1913
    y2i := x1i + x2i;
1914
    y6r := x1r - x2r;
1915
    y6i := x1i - x2i;
1916
    x0r := a[ 4 ] + a[ 21 ];
1917
    x0i := a[ 5 ] - a[ 20 ];
1918
    x1r := wk2i * x0r - wk2r * x0i;
1919
    x1i := wk2i * x0i + wk2r * x0r;
1920
    x0r := a[ 12 ] + a[ 29 ];
1921
    x0i := a[ 13 ] - a[ 28 ];
1922
    x2r := wk2r * x0r - wk2i * x0i;
1923
    x2i := wk2r * x0i + wk2i * x0r;
1924
    y10r := x1r - x2r;
1925
    y10i := x1i - x2i;
1926
    y14r := x1r + x2r;
1927
    y14i := x1i + x2i;
1928
    x0r := a[ 6 ] - a[ 23 ];
1929
    x0i := a[ 7 ] + a[ 22 ];
1930
    x1r := wk3r * x0r - wk3i * x0i;
1931
    x1i := wk3r * x0i + wk3i * x0r;
1932
    x0r := a[ 14 ] - a[ 31 ];
1933
    x0i := a[ 15 ] + a[ 30 ];
1934
    x2r := wk1i * x0r - wk1r * x0i;
1935
    x2i := wk1i * x0i + wk1r * x0r;
1936
    y3r := x1r + x2r;
1937
    y3i := x1i + x2i;
1938
    y7r := x1r - x2r;
1939
    y7i := x1i - x2i;
1940
    x0r := a[ 6 ] + a[ 23 ];
1941
    x0i := a[ 7 ] - a[ 22 ];
1942
    x1r := wk1i * x0r + wk1r * x0i;
1943
    x1i := wk1i * x0i - wk1r * x0r;
1944
    x0r := a[ 14 ] + a[ 31 ];
1945
    x0i := a[ 15 ] - a[ 30 ];
1946
    x2r := wk3i * x0r - wk3r * x0i;
1947
    x2i := wk3i * x0i + wk3r * x0r;
1948
    y11r := x1r + x2r;
1949
    y11i := x1i + x2i;
1950
    y15r := x1r - x2r;
1951
    y15i := x1i - x2i;
1952
    x1r := y0r + y2r;
1953
    x1i := y0i + y2i;
1954
    x2r := y1r + y3r;
1955
    x2i := y1i + y3i;
1956
    a[ 0 ] := x1r + x2r;
1957
    a[ 1 ] := x1i + x2i;
1958
    a[ 2 ] := x1r - x2r;
1959
    a[ 3 ] := x1i - x2i;
1960
    x1r := y0r - y2r;
1961
    x1i := y0i - y2i;
1962
    x2r := y1r - y3r;
1963
    x2i := y1i - y3i;
1964
    a[ 4 ] := x1r - x2i;
1965
    a[ 5 ] := x1i + x2r;
1966
    a[ 6 ] := x1r + x2i;
1967
    a[ 7 ] := x1i - x2r;
1968
    x1r := y4r - y6i;
1969
    x1i := y4i + y6r;
1970
    x0r := y5r - y7i;
1971
    x0i := y5i + y7r;
1972
    x2r := wn4r * ( x0r - x0i );
1973
    x2i := wn4r * ( x0i + x0r );
1974
    a[ 8 ] := x1r + x2r;
1975
    a[ 9 ] := x1i + x2i;
1976
    a[ 10 ] := x1r - x2r;
1977
    a[ 11 ] := x1i - x2i;
1978
    x1r := y4r + y6i;
1979
    x1i := y4i - y6r;
1980
    x0r := y5r + y7i;
1981
    x0i := y5i - y7r;
1982
    x2r := wn4r * ( x0r - x0i );
1983
    x2i := wn4r * ( x0i + x0r );
1984
    a[ 12 ] := x1r - x2i;
1985
    a[ 13 ] := x1i + x2r;
1986
    a[ 14 ] := x1r + x2i;
1987
    a[ 15 ] := x1i - x2r;
1988
    x1r := y8r + y10r;
1989
    x1i := y8i + y10i;
1990
    x2r := y9r - y11r;
1991
    x2i := y9i - y11i;
1992
    a[ 16 ] := x1r + x2r;
1993
    a[ 17 ] := x1i + x2i;
1994
    a[ 18 ] := x1r - x2r;
1995
    a[ 19 ] := x1i - x2i;
1996
    x1r := y8r - y10r;
1997
    x1i := y8i - y10i;
1998
    x2r := y9r + y11r;
1999
    x2i := y9i + y11i;
2000
    a[ 20 ] := x1r - x2i;
2001
    a[ 21 ] := x1i + x2r;
2002
    a[ 22 ] := x1r + x2i;
2003
    a[ 23 ] := x1i - x2r;
2004
    x1r := y12r - y14i;
2005
    x1i := y12i + y14r;
2006
    x0r := y13r + y15i;
2007
    x0i := y13i - y15r;
2008
    x2r := wn4r * ( x0r - x0i );
2009
    x2i := wn4r * ( x0i + x0r );
2010
    a[ 24 ] := x1r + x2r;
2011
    a[ 25 ] := x1i + x2i;
2012
    a[ 26 ] := x1r - x2r;
2013
    a[ 27 ] := x1i - x2i;
2014
    x1r := y12r + y14i;
2015
    x1i := y12i - y14r;
2016
    x0r := y13r - y15i;
2017
    x0i := y13i + y15r;
2018
    x2r := wn4r * ( x0r - x0i );
2019
    x2i := wn4r * ( x0i + x0r );
2020
    a[ 28 ] := x1r - x2i;
2021
    a[ 29 ] := x1i + x2r;
2022
    a[ 30 ] := x1r + x2i;
2023
    a[ 31 ] := x1i - x2r;
2024
end;
2025

2026
////////////////////////////////////////////////////////////////////////////////////////////////////
2027

2028
procedure TDiscreteTrans.cftf081( var a:array of Double; var w:array of Double );
2029
var
2030
   wn4r, x0r, x0i, x1r, x1i, x2r, x2i, x3r, x3i,
2031
   y0r, y0i, y1r, y1i, y2r, y2i, y3r, y3i,
2032
   y4r, y4i, y5r, y5i, y6r, y6i, y7r, y7i :Double;
2033
begin
2034
     wn4r := w[ 1 ];
2035

2036
     x0r := a[ 00 ] + a[ 08 ];    x0i := a[ 01 ] + a[ 09 ];
2037
     x1r := a[ 00 ] - a[ 08 ];    x1i := a[ 01 ] - a[ 09 ];
2038
     x2r := a[ 04 ] + a[ 12 ];    x2i := a[ 05 ] + a[ 13 ];
2039
     x3r := a[ 04 ] - a[ 12 ];    x3i := a[ 05 ] - a[ 13 ];
2040
     y0r := x0r + x2r;            y0i := x0i + x2i;
2041
     y2r := x0r - x2r;            y2i := x0i - x2i;
2042
     y1r := x1r - x3i;            y1i := x1i + x3r;
2043
     y3r := x1r + x3i;            y3i := x1i - x3r;
2044
     x0r := a[ 02 ] + a[ 10 ];    x0i := a[ 03 ] + a[ 11 ];
2045
     x1r := a[ 02 ] - a[ 10 ];    x1i := a[ 03 ] - a[ 11 ];
2046
     x2r := a[ 06 ] + a[ 14 ];    x2i := a[ 07 ] + a[ 15 ];
2047
     x3r := a[ 06 ] - a[ 14 ];    x3i := a[ 07 ] - a[ 15 ];
2048
     y4r := x0r + x2r;            y4i := x0i + x2i;
2049
     y6r := x0r - x2r;            y6i := x0i - x2i;
2050
     x0r := x1r - x3i;            x0i := x1i + x3r;
2051
     x2r := x1r + x3i;            x2i := x1i - x3r;
2052
     y5r := wn4r * ( x0r - x0i ); y5i := wn4r * ( x0r + x0i );
2053
     y7r := wn4r * ( x2r - x2i ); y7i := wn4r * ( x2r + x2i );
2054

2055
     a[ 08 ] := y1r + y5r;
2056
     a[ 09 ] := y1i + y5i;
2057
     a[ 10 ] := y1r - y5r;
2058
     a[ 11 ] := y1i - y5i;
2059
     a[ 12 ] := y3r - y7i;
2060
     a[ 13 ] := y3i + y7r;
2061
     a[ 14 ] := y3r + y7i;
2062
     a[ 15 ] := y3i - y7r;
2063
     a[ 00 ] := y0r + y4r;
2064
     a[ 01 ] := y0i + y4i;
2065
     a[ 02 ] := y0r - y4r;
2066
     a[ 03 ] := y0i - y4i;
2067
     a[ 04 ] := y2r - y6i;
2068
     a[ 05 ] := y2i + y6r;
2069
     a[ 06 ] := y2r + y6i;
2070
     a[ 07 ] := y2i - y6r;
2071
end;
2072

2073
procedure TDiscreteTrans.cftf082( var a:array of Double; var w:array of Double );
2074
var
2075
   wn4r, wk1r, wk1i, x0r, x0i, x1r, x1i,
2076
   y0r, y0i, y1r, y1i, y2r, y2i, y3r, y3i,
2077
   y4r, y4i, y5r, y5i, y6r, y6i, y7r, y7i :Double;
2078
begin
2079
    wn4r := w[ 1 ];
2080
    wk1r := w[ 2 ];
2081
    wk1i := w[ 3 ];
2082

2083
    y0r := a[ 00 ] - a[ 09 ];       y0i := a[ 01 ] + a[ 08 ];
2084
    y1r := a[ 00 ] + a[ 09 ];       y1i := a[ 01 ] - a[ 08 ];
2085
    x0r := a[ 04 ] - a[ 13 ];       x0i := a[ 05 ] + a[ 12 ];
2086
    y2r := wn4r * ( x0r - x0i );    y2i := wn4r * ( x0i + x0r );
2087
    x0r := a[ 04 ] + a[ 13 ];       x0i := a[ 05 ] - a[ 12 ];
2088
    y3r := wn4r * ( x0r - x0i );    y3i := wn4r * ( x0i + x0r );
2089
    x0r := a[ 02 ] - a[ 11 ];       x0i := a[ 03 ] + a[ 10 ];
2090
    y4r := wk1r * x0r - wk1i * x0i; y4i := wk1r * x0i + wk1i * x0r;
2091
    x0r := a[ 02 ] + a[ 11 ];       x0i := a[ 03 ] - a[ 10 ];
2092
    y5r := wk1i * x0r - wk1r * x0i; y5i := wk1i * x0i + wk1r * x0r;
2093
    x0r := a[ 06 ] - a[ 15 ];       x0i := a[ 07 ] + a[ 14 ];
2094
    y6r := wk1i * x0r - wk1r * x0i; y6i := wk1i * x0i + wk1r * x0r;
2095
    x0r := a[ 06 ] + a[ 15 ];       x0i := a[ 07 ] - a[ 14 ];
2096
    y7r := wk1r * x0r - wk1i * x0i; y7i := wk1r * x0i + wk1i * x0r;
2097
    x0r := y0r + y2r;               x0i := y0i + y2i;
2098
    x1r := y4r + y6r;               x1i := y4i + y6i;
2099
    a[ 00 ] := x0r + x1r;
2100
    a[ 01 ] := x0i + x1i;
2101
    a[ 02 ] := x0r - x1r;
2102
    a[ 03 ] := x0i - x1i;
2103

2104
    x0r := y0r - y2r; x0i := y0i - y2i;
2105
    x1r := y4r - y6r; x1i := y4i - y6i;
2106
    a[ 04 ] := x0r - x1i;
2107
    a[ 05 ] := x0i + x1r;
2108
    a[ 06 ] := x0r + x1i;
2109
    a[ 07 ] := x0i - x1r;
2110

2111
    x0r := y1r - y3i; x0i := y1i + y3r;
2112
    x1r := y5r - y7r; x1i := y5i - y7i;
2113
    a[ 08 ] := x0r + x1r;
2114
    a[ 09 ] := x0i + x1i;
2115
    a[ 10 ] := x0r - x1r;
2116
    a[ 11 ] := x0i - x1i;
2117

2118
    x0r := y1r + y3i; x0i := y1i - y3r;
2119
    x1r := y5r + y7r; x1i := y5i + y7i;
2120
    a[ 12 ] := x0r - x1i;
2121
    a[ 13 ] := x0i + x1r;
2122
    a[ 14 ] := x0r + x1i;
2123
    a[ 15 ] := x0i - x1r;
2124
end;
2125

2126
////////////////////////////////////////////////////////////////////////////////////////////////////
2127

2128
procedure TDiscreteTrans.cftleaf( const n,isplt:Integer; var a:array of Double );
2129
begin
2130
     if n = 512 then
2131
     begin
2132
          cftmdl1( 128, a       , _W[ _NW -  64 ] );
2133
          cftf161(      a       , _W[ _NW -   8 ] );
2134
          cftf162(      a[  32 ], _W[ _NW -  32 ] );
2135
          cftf161(      a[  64 ], _W[ _NW -   8 ] );
2136
          cftf161(      a[  96 ], _W[ _NW -   8 ] );
2137
          cftmdl2( 128, a[ 128 ], _W[ _NW - 128 ] );
2138
          cftf161(      a[ 128 ], _W[ _NW -   8 ] );
2139
          cftf162(      a[ 160 ], _W[ _NW -  32 ] );
2140
          cftf161(      a[ 192 ], _W[ _NW -   8 ] );
2141
          cftf162(      a[ 224 ], _W[ _NW -  32 ] );
2142
          cftmdl1( 128, a[ 256 ], _W[ _NW -  64 ] );
2143
          cftf161(      a[ 256 ], _W[ _NW -   8 ] );
2144
          cftf162(      a[ 288 ], _W[ _NW -  32 ] );
2145
          cftf161(      a[ 320 ], _W[ _NW -   8 ] );
2146
          cftf161(      a[ 352 ], _W[ _NW -   8 ] );
2147

2148
          if isplt <> 0 then
2149
          begin
2150
               cftmdl1( 128, a[ 384 ], _W[ _NW - 64 ] );
2151
               cftf161(      a[ 480 ], _W[ _NW -  8 ] );
2152
          end
2153
          else
2154
          begin
2155
               cftmdl2( 128, a[ 384 ], _W[ _NW - 128 ] );
2156
               cftf162(      a[ 480 ], _W[ _NW -  32 ] );
2157
          end;
2158

2159
          cftf161( a[ 384 ], _W[ _NW -  8 ] );
2160
          cftf162( a[ 416 ], _W[ _NW - 32 ] );
2161
          cftf161( a[ 448 ], _W[ _NW -  8 ] );
2162
     end
2163
     else
2164
     begin
2165
          cftmdl1( 64, a       , _W[ _NW - 32 ] );
2166
          cftf081(     a       , _W[ _NW -  8 ] );
2167
          cftf082(     a[  16 ], _W[ _NW -  8 ] );
2168
          cftf081(     a[  32 ], _W[ _NW -  8 ] );
2169
          cftf081(     a[  48 ], _W[ _NW -  8 ] );
2170
          cftmdl2( 64, a[  64 ], _W[ _NW - 64 ] );
2171
          cftf081(     a[  64 ], _W[ _NW -  8 ] );
2172
          cftf082(     a[  80 ], _W[ _NW -  8 ] );
2173
          cftf081(     a[  96 ], _W[ _NW -  8 ] );
2174
          cftf082(     a[ 112 ], _W[ _NW -  8 ] );
2175
          cftmdl1( 64, a[ 128 ], _W[ _NW - 32 ] );
2176
          cftf081(     a[ 128 ], _W[ _NW -  8 ] );
2177
          cftf082(     a[ 144 ], _W[ _NW -  8 ] );
2178
          cftf081(     a[ 160 ], _W[ _NW -  8 ] );
2179
          cftf081(     a[ 176 ], _W[ _NW -  8 ] );
2180

2181
          if isplt <> 0 then
2182
          begin
2183
               cftmdl1( 64, a[ 192 ], _W[ _NW - 32 ] );
2184
               cftf081(     a[ 240 ], _W[ _NW -  8 ] );
2185
          end
2186
          else
2187
          begin
2188
               cftmdl2( 64, a[ 192 ], _W[ _NW - 64 ] );
2189
               cftf082(     a[ 240 ], _W[ _NW -  8 ] );
2190
          end;
2191

2192
          cftf081( a[ 192 ], _W[ _NW - 8 ] );
2193
          cftf082( a[ 208 ], _W[ _NW - 8 ] );
2194
          cftf081( a[ 224 ], _W[ _NW - 8 ] );
2195
    end
2196
end;
2197

2198
////////////////////////////////////////////////////////////////////////////////////////////////////
2199

2200
procedure TDiscreteTrans.cftrec4;
2201
var
2202
   isplt, j, k, m :Integer;
2203
begin
2204
     m := _TempN;
2205
     while m > 512 do
2206
     begin
2207
          m := m shr 2;
2208
          cftmdl1( m, _TEMP[ _TempN - m ], _W[ _NW - ( m shr 1 ) ])
2209
     end;
2210

2211
     cftleaf( m, 1, _TEMP[ _TempN - m ] );
2212

2213
     k := 0;
2214
     j := _TempN - m;
2215
     while j > 0 do
2216
     begin
2217
          Inc( k );
2218

2219
          isplt := cfttree( m, j, k );
2220
          cftleaf( m, isplt, _TEMP[ j - m ] );
2221

2222
          Dec( j, m )
2223
     end
2224
end;
2225

2226
////////////////////////////////////////////////////////////////////////////////////////////////////
2227

2228
procedure TDiscreteTrans.cftfx41;
2229
begin
2230
    if _TempN = 128 then
2231
    begin
2232
         cftf161( _TEMP      , _W[ _NW -  8 ] );
2233
         cftf162( _TEMP[ 32 ], _W[ _NW - 32 ] );
2234
         cftf161( _TEMP[ 64 ], _W[ _NW -  8 ] );
2235
         cftf161( _TEMP[ 96 ], _W[ _NW -  8 ] );
2236
    end
2237
    else
2238
    begin
2239
         cftf081( _TEMP      , _W[ _NW -  8 ] );
2240
         cftf082( _TEMP[ 16 ], _W[ _NW -  8 ] );
2241
         cftf081( _TEMP[ 32 ], _W[ _NW -  8 ] );
2242
         cftf081( _TEMP[ 48 ], _W[ _NW -  8 ] );
2243
    end
2244
end;
2245

2246
////////////////////////////////////////////////////////////////////////////////////////////////////
2247

2248
procedure TDiscreteTrans.cftf040;
2249
var
2250
   x0r, x0i,
2251
   x1r, x1i,
2252
   x2r, x2i,
2253
   x3r, x3i :Double;
2254
begin
2255
     x0r := _TEMP[ 0 ] + _TEMP[ 4 ];  x0i := _TEMP[ 1 ] + _TEMP[ 5 ];
2256
     x1r := _TEMP[ 0 ] - _TEMP[ 4 ];  x1i := _TEMP[ 1 ] - _TEMP[ 5 ];
2257
     x2r := _TEMP[ 2 ] + _TEMP[ 6 ];  x2i := _TEMP[ 3 ] + _TEMP[ 7 ];
2258
     x3r := _TEMP[ 2 ] - _TEMP[ 6 ];  x3i := _TEMP[ 3 ] - _TEMP[ 7 ];
2259

2260
     _TEMP[ 0 ] := x0r + x2r;  _TEMP[ 1 ] := x0i + x2i;
2261
     _TEMP[ 2 ] := x1r - x3i;  _TEMP[ 3 ] := x1i + x3r;
2262
     _TEMP[ 4 ] := x0r - x2r;  _TEMP[ 5 ] := x0i - x2i;
2263
     _TEMP[ 6 ] := x1r + x3i;  _TEMP[ 7 ] := x1i - x3r;
2264
end;
2265

2266
procedure TDiscreteTrans.cftb040;
2267
var
2268
   x0r, x0i,
2269
   x1r, x1i,
2270
   x2r, x2i,
2271
   x3r, x3i :Double;
2272
begin
2273
     x0r := _TEMP[ 0 ] + _TEMP[ 4 ];  x0i := _TEMP[ 1 ] + _TEMP[ 5 ];
2274
     x1r := _TEMP[ 0 ] - _TEMP[ 4 ];  x1i := _TEMP[ 1 ] - _TEMP[ 5 ];
2275
     x2r := _TEMP[ 2 ] + _TEMP[ 6 ];  x2i := _TEMP[ 3 ] + _TEMP[ 7 ];
2276
     x3r := _TEMP[ 2 ] - _TEMP[ 6 ];  x3i := _TEMP[ 3 ] - _TEMP[ 7 ];
2277

2278
     _TEMP[ 0 ] := x0r + x2r;  _TEMP[ 1 ] := x0i + x2i;
2279
     _TEMP[ 2 ] := x1r + x3i;  _TEMP[ 3 ] := x1i - x3r;
2280
     _TEMP[ 4 ] := x0r - x2r;  _TEMP[ 5 ] := x0i - x2i;
2281
     _TEMP[ 6 ] := x1r - x3i;  _TEMP[ 7 ] := x1i + x3r;
2282
end;
2283

2284
////////////////////////////////////////////////////////////////////////////////////////////////////
2285

2286
procedure TDiscreteTrans.cftx020;
2287
var
2288
   x0r, x0i :Double;
2289
begin
2290
            x0r := _TEMP[ 0 ] - _TEMP[ 2 ];
2291
            x0i := _TEMP[ 1 ] - _TEMP[ 3 ];
2292

2293
     _TEMP[ 0 ] := _TEMP[ 0 ] + _TEMP[ 2 ];
2294
     _TEMP[ 1 ] := _TEMP[ 1 ] + _TEMP[ 3 ];
2295

2296
     _TEMP[ 2 ] := x0r;
2297
     _TEMP[ 3 ] := x0i;
2298
end;
2299

2300
//&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& protected
2301

2302
//XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX アクセス
2303

2304
procedure TDiscreteTrans.SetCount( const Count_:Integer );
2305
begin
2306
     _Count := Count_;
2307
end;
2308

2309
////////////////////////////////////////////////////////////////////////////////////////////////////
2310

2311
procedure TDiscreteTrans.cftfsub;
2312
begin
2313
     if _TempN > 8 then
2314
     begin
2315
          if _TempN > 32 then
2316
          begin
2317
               cftf1st( _W[ _NW - ( _TempN shr 2  ) ] );
2318

2319
               if _TempN > 512 then cftrec4
2320
                               else
2321
               if _TempN > 128 then cftleaf( _TempN, 1, _TEMP )
2322
                               else cftfx41;
2323

2324
               bitrv2;
2325
          end
2326
          else
2327
          if _TempN = 32 then
2328
          begin
2329
               cftf161( _TEMP, _W[ _NW - 8 ] );
2330
               bitrv216;
2331
          end
2332
          else
2333
          begin
2334
               cftf081( _TEMP, _W );
2335
               bitrv208;
2336
          end
2337
     end
2338
     else
2339
     if _TempN = 8 then cftf040
2340
     else
2341
     if _TempN = 4 then cftx020
2342
end;
2343

2344
procedure TDiscreteTrans.cftbsub;
2345
begin
2346
     if _TempN > 8 then
2347
     begin
2348
          if _TempN > 32 then
2349
          begin
2350
               cftb1st( _W[ _NW - ( _TempN shr 2  ) ] );
2351

2352
               if _TempN > 512 then cftrec4
2353
                               else
2354
               if _TempN > 128 then cftleaf( _TempN, 1, _TEMP )
2355
                               else cftfx41;
2356

2357
               bitrv2conj;
2358
          end
2359
          else
2360
          if _TempN = 32 then
2361
          begin
2362
               cftf161( _TEMP, _W[ _NW - 8 ] );
2363

2364
               bitrv216neg;
2365
          end
2366
          else
2367
          begin
2368
               cftf081( _TEMP, _W );
2369

2370
               bitrv208neg;
2371
          end
2372
     end
2373
     else
2374
     if _TempN = 8 then cftb040
2375
     else
2376
     if _TempN = 4 then cftx020
2377
end;
2378

2379
////////////////////////////////////////////////////////////////////////////////////////////////////
2380

2381
procedure TDiscreteTrans.rftfsub( var c:array of Double );
2382
var
2383
   j, k, kk, ks, m :Integer;
2384
   wkr, wki,
2385
    xr,  xi,
2386
    yr,  yi :Double;
2387
begin
2388
     m := _TempN shr 1;
2389

2390
     ks := 2 * _NC div m;
2391
     kk := 0;
2392
     j := 2;
2393
     while j < m do
2394
     begin
2395
          k := _TempN - j;
2396
          Inc( kk, ks );
2397

2398
          wkr := 0.5 - c[ _NC - kk ];
2399
          wki :=       c[       kk ];
2400

2401
          xr := _TEMP[ j     ] - _TEMP[ k     ];
2402
          xi := _TEMP[ j + 1 ] + _TEMP[ k + 1 ];
2403

2404
          yr := wkr * xr - wki * xi;
2405
          yi := wkr * xi + wki * xr;
2406

2407
          _TEMP[ j     ] := _TEMP[ j     ] - yr;
2408
          _TEMP[ j + 1 ] := _TEMP[ j + 1 ] - yi;
2409
          _TEMP[ k     ] := _TEMP[ k     ] + yr;
2410
          _TEMP[ k + 1 ] := _TEMP[ k + 1 ] - yi;
2411

2412
          Inc( j, 2 )
2413
     end
2414
end;
2415

2416
procedure TDiscreteTrans.rftbsub( var c:array of Double );
2417
var
2418
   j, k, kk, ks, m :Integer;
2419
   wkr, wki,
2420
    xr,  xi,
2421
    yr,  yi :Double;
2422
begin
2423
     m := _TempN shr 1;
2424

2425
     ks := 2 * _NC div m;
2426
     kk := 0;
2427
     j := 2;
2428
     while j < m do
2429
     begin
2430
          k  := _TempN - j;
2431
          kk := kk + ks;
2432

2433
          wkr := 0.5 - c[ _NC - kk ];
2434
          wki :=       c[       kk ];
2435

2436
          xr := _TEMP[ j     ] - _TEMP[ k     ];
2437
          xi := _TEMP[ j + 1 ] + _TEMP[ k + 1 ];
2438

2439
          yr := wkr * xr + wki * xi;
2440
          yi := wkr * xi - wki * xr;
2441

2442
          _TEMP[ j     ] := _TEMP[ j     ] - yr;
2443
          _TEMP[ j + 1 ] := _TEMP[ j + 1 ] - yi;
2444
          _TEMP[ k     ] := _TEMP[ k     ] + yr;
2445
          _TEMP[ k + 1 ] := _TEMP[ k + 1 ] - yi;
2446

2447
          Inc( j, 2 )
2448
     end
2449
end;
2450

2451
////////////////////////////////////////////////////////////////////////////////////////////////////
2452

2453
procedure TDiscreteTrans.Normalize;
2454
var
2455
   N :Integer;
2456
   P :PDouble;
2457
begin
2458
     P := @_TEMP[ 0 ];
2459
     for N := 1 to _TempN do
2460
     begin
2461
          P^ := _NormW * P^; Inc( P )
2462
     end
2463
end;
2464

2465

2466
////////////////////////////////////////////////////////////////////////////////////////////////////
2467

2468
procedure TDiscreteTrans.MakeTableW;
2469
begin
2470
     _NW := _IP[ 0 ];
2471

2472
     if _TempN > _NW shl 2 then
2473
     begin
2474
          _NW := _TempN shr 2;
2475

2476
          makewt
2477
     end
2478
end;
2479

2480
procedure TDiscreteTrans.MakeTableC;
2481
begin
2482
     _NC := _IP[ 1 ];
2483

2484
     if _TempN > _NC shl 2 then
2485
     begin
2486
          _NC := _TempN shr 2;
2487

2488
          makect( _W[ _NW ] )
2489
     end
2490
end;
2491

2492
//&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& public
2493

2494
constructor TDiscreteTrans.Create;
2495
begin
2496
     inherited;
2497

2498
end;
2499

2500
destructor TDiscreteTrans.Destroy;
2501
begin
2502

2503
     inherited;
2504
end;
2505

2506

2507
//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$【ルーチン】
2508

2509
//################################################################################################## □
2510

2511
initialization //############################################################################ 初期化
2512

2513
finalization //############################################################################## 終了化
2514

2515
end. //############################################################################################# ■
2516