1 /* temperton/gpfa5f.f -- translated by f2c (version 20050501).
2 You must link the resulting object file with libf2c:
3 on Microsoft Windows system, link with libf2c.lib;
4 on Linux or Unix systems, link with .../path/to/libf2c.a -lm
5 or, if you install libf2c.a in a standard place, with -lf2c -lm
6 -- in that order, at the end of the command line, as in
7 cc *.o -lf2c -lm
8 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
9
10 http://www.netlib.org/f2c/libf2c.zip
11 */
12
13 #ifdef __cplusplus
14 extern "C" {
15 #endif
16 #include "v3p_netlib.h"
17
18 /* Table of constant values */
19
20 static integer c__5 = 5;
21
22 /* fortran version of *gpfa5* - */
23 /* radix-5 section of self-sorting, in-place, */
24 /* generalized pfa */
25
26 /* ------------------------------------------------------------------- */
27
28 /*< subroutine gpfa5f(a,b,trigs,inc,jump,n,mm,lot,isign) >*/
gpfa5f_(real * a,real * b,real * trigs,integer * inc,integer * jump,integer * n,integer * mm,integer * lot,integer * isign)29 /* Subroutine */ int gpfa5f_(real *a, real *b, real *trigs, integer *inc,
30 integer *jump, integer *n, integer *mm, integer *lot, integer *isign)
31 {
32 /* Initialized data */
33
34 static real sin36 = (float).587785252292473; /* constant */
35 static real sin72 = (float).951056516295154; /* constant */
36 static real qrt5 = (float).559016994374947; /* constant */
37 static integer lvr = 128; /* constant */
38
39 /* System generated locals */
40 integer i__1, i__2, i__3, i__4, i__5, i__6, i__7, i__8, i__9, i__10;
41
42 /* Builtin functions */
43 integer pow_ii(integer *, integer *);
44
45 /* Local variables */
46 integer j, k, l, m;
47 real s, c1, c2, c3;
48 integer n5;
49 real t1, t2, t3, t4, t5, t6, t7, t8, t9, u1, u2, u3, u4, u5, u6, u7, u8,
50 u9;
51 integer ja, jb, la, jc, jd, nb, je, jf, jg, jh;
52 real t10, t11, u10, u11, ax, bx;
53 integer mh, kk, ll, ji, jj, jk, mu, nu, jl, jm, jn, jo, jp, jq, jr, js,
54 jt, ju, jv, jw, jx, jy;
55 real co1=0, co2=0, co3=0, co4=0, si1=0, si2=0, si3=0, si4=0,
56 aja, ajb, ajc, ajd, aje, bjb,
57 bje, bjc, bjd, bja, ajf, ajk, bjf, bjk, ajg, ajj, ajh, aji, ajl,
58 ajq, bjg, bjj, bjh, bji, bjl, bjq, ajo, ajm, ajn, ajr, ajw, bjo,
59 bjm, bjn, bjr, bjw, ajt, ajs, ajx, ajp, bjt, bjs, bjx, bjp, ajv,
60 ajy, aju, bjv, bjy, bju;
61 integer inq, ink, jjj, ninc, left, nvex, ipass, nblox, jstep, laincl,
62 jstepl, istart, jstepx;
63
64 /*< real a(*), b(*), trigs(*) >*/
65 /*< integer inc, jump, n, mm, lot, isign >*/
66 /*< real s, ax, bx, c1, c2, c3 >*/
67 /*< real t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11 >*/
68 /*< real u1, u2, u3, u4, u5, u6, u7, u8, u9, u10, u11 >*/
69 /*< real co1, co2, co3, co4, si1, si2, si3, si4 >*/
70 /*< real aja, ajb, ajc, ajd, aje, bjb, bje, bjc >*/
71 /*< real bjd, bja, ajf, ajk, bjf, bjk, ajg, ajj >*/
72 /*< real ajh, aji, ajl, ajq, bjg, bjj, bjh, bji >*/
73 /*< real bjl, bjq, ajo, ajm, ajn, ajr, ajw, bjo >*/
74 /*< real bjm, bjn, bjr, bjw, ajt, ajs, ajx, ajp >*/
75 /*< real bjt, bjs, bjx, bjp, ajv, ajy, aju, bjv >*/
76 /*< real bjy, bju >*/
77 /*< >*/
78 /* Parameter adjustments */
79 --trigs;
80 --b;
81 --a;
82
83 /* Function Body */
84 /*< data lvr/128/ >*/
85
86 /* *************************************************************** */
87 /* * * */
88 /* * N.B. LVR = LENGTH OF VECTOR REGISTERS, SET TO 128 FOR C90. * */
89 /* * RESET TO 64 FOR OTHER CRAY MACHINES, OR TO ANY LARGE VALUE * */
90 /* * (GREATER THAN OR EQUAL TO LOT) FOR A SCALAR COMPUTER. * */
91 /* * * */
92 /* *************************************************************** */
93
94 /*< n5 = 5 ** mm >*/
95 n5 = pow_ii(&c__5, mm);
96 /*< inq = n / n5 >*/
97 inq = *n / n5;
98 /*< jstepx = (n5-n) * inc >*/
99 jstepx = (n5 - *n) * *inc;
100 /*< ninc = n * inc >*/
101 ninc = *n * *inc;
102 /*< ink = inc * inq >*/
103 ink = *inc * inq;
104 /*< mu = mod(inq,5) >*/
105 mu = inq % 5;
106 /*< if (isign.eq.-1) mu = 5 - mu >*/
107 if (*isign == -1) {
108 mu = 5 - mu;
109 }
110
111 /*< m = mm >*/
112 m = *mm;
113 /*< mh = (m+1)/2 >*/
114 mh = (m + 1) / 2;
115 /*< s = float(isign) >*/
116 // s = (real) (*isign);
117 /*< c1 = qrt5 >*/
118 c1 = qrt5;
119 /*< c2 = sin72 >*/
120 c2 = sin72;
121 /*< c3 = sin36 >*/
122 c3 = sin36;
123 /*< if (mu.eq.2.or.mu.eq.3) then >*/
124 if (mu == 2 || mu == 3) {
125 /*< c1 = -c1 >*/
126 c1 = -c1;
127 /*< c2 = sin36 >*/
128 c2 = sin36;
129 /*< c3 = sin72 >*/
130 c3 = sin72;
131 /*< endif >*/
132 }
133 /*< if (mu.eq.3.or.mu.eq.4) c2 = -c2 >*/
134 if (mu == 3 || mu == 4) {
135 c2 = -c2;
136 }
137 /*< if (mu.eq.2.or.mu.eq.4) c3 = -c3 >*/
138 if (mu == 2 || mu == 4) {
139 c3 = -c3;
140 }
141
142 /*< nblox = 1 + (lot-1)/lvr >*/
143 nblox = (*lot - 1) / lvr + 1;
144 /*< left = lot >*/
145 left = *lot;
146 /*< s = float(isign) >*/
147 s = (real) (*isign);
148 /*< istart = 1 >*/
149 istart = 1;
150
151 /* loop on blocks of lvr transforms */
152 /* -------------------------------- */
153 /*< do 500 nb = 1 , nblox >*/
154 i__1 = nblox;
155 for (nb = 1; nb <= i__1; ++nb) {
156
157 /*< if (left.le.lvr) then >*/
158 if (left <= lvr) {
159 /*< nvex = left >*/
160 nvex = left;
161 /*< else if (left.lt.(2*lvr)) then >*/
162 } else if (left < lvr << 1) {
163 /*< nvex = left/2 >*/
164 nvex = left / 2;
165 /*< nvex = nvex + mod(nvex,2) >*/
166 nvex += nvex % 2;
167 /*< else >*/
168 } else {
169 /*< nvex = lvr >*/
170 nvex = lvr;
171 /*< endif >*/
172 }
173 /*< left = left - nvex >*/
174 left -= nvex;
175
176 /*< la = 1 >*/
177 la = 1;
178
179 /* loop on type I radix-5 passes */
180 /* ----------------------------- */
181 /*< do 160 ipass = 1 , mh >*/
182 i__2 = mh;
183 for (ipass = 1; ipass <= i__2; ++ipass) {
184 /*< jstep = (n*inc) / (5*la) >*/
185 jstep = *n * *inc / (la * 5);
186 /*< jstepl = jstep - ninc >*/
187 jstepl = jstep - ninc;
188 /*< kk = 0 >*/
189 kk = 0;
190
191 /* loop on k */
192 /* --------- */
193 /*< do 150 k = 0 , jstep-ink , ink >*/
194 i__3 = jstep - ink;
195 i__4 = ink;
196 for (k = 0; i__4 < 0 ? k >= i__3 : k <= i__3; k += i__4) {
197
198 /*< if (k.gt.0) then >*/
199 if (k > 0) {
200 /*< co1 = trigs(kk+1) >*/
201 co1 = trigs[kk + 1];
202 /*< si1 = s*trigs(kk+2) >*/
203 si1 = s * trigs[kk + 2];
204 /*< co2 = trigs(2*kk+1) >*/
205 co2 = trigs[(kk << 1) + 1];
206 /*< si2 = s*trigs(2*kk+2) >*/
207 si2 = s * trigs[(kk << 1) + 2];
208 /*< co3 = trigs(3*kk+1) >*/
209 co3 = trigs[kk * 3 + 1];
210 /*< si3 = s*trigs(3*kk+2) >*/
211 si3 = s * trigs[kk * 3 + 2];
212 /*< co4 = trigs(4*kk+1) >*/
213 co4 = trigs[(kk << 2) + 1];
214 /*< si4 = s*trigs(4*kk+2) >*/
215 si4 = s * trigs[(kk << 2) + 2];
216 /*< endif >*/
217 }
218
219 /* loop along transform */
220 /* -------------------- */
221 /*< do 140 jjj = k , (n-1)*inc , 5*jstep >*/
222 i__5 = (*n - 1) * *inc;
223 i__6 = jstep * 5;
224 for (jjj = k; i__6 < 0 ? jjj >= i__5 : jjj <= i__5; jjj +=
225 i__6) {
226 /*< ja = istart + jjj >*/
227 ja = istart + jjj;
228
229 /* "transverse" loop */
230 /* ----------------- */
231 /*< do 135 nu = 1 , inq >*/
232 i__7 = inq;
233 for (nu = 1; nu <= i__7; ++nu) {
234 /*< jb = ja + jstepl >*/
235 jb = ja + jstepl;
236 /*< if (jb.lt.istart) jb = jb + ninc >*/
237 if (jb < istart) {
238 jb += ninc;
239 }
240 /*< jc = jb + jstepl >*/
241 jc = jb + jstepl;
242 /*< if (jc.lt.istart) jc = jc + ninc >*/
243 if (jc < istart) {
244 jc += ninc;
245 }
246 /*< jd = jc + jstepl >*/
247 jd = jc + jstepl;
248 /*< if (jd.lt.istart) jd = jd + ninc >*/
249 if (jd < istart) {
250 jd += ninc;
251 }
252 /*< je = jd + jstepl >*/
253 je = jd + jstepl;
254 /*< if (je.lt.istart) je = je + ninc >*/
255 if (je < istart) {
256 je += ninc;
257 }
258 /*< j = 0 >*/
259 j = 0;
260
261 /* loop across transforms */
262 /* ---------------------- */
263 /*< if (k.eq.0) then >*/
264 if (k == 0) {
265
266 /* dir$ ivdep, shortloop */
267 /*< do 110 l = 1 , nvex >*/
268 i__8 = nvex;
269 for (l = 1; l <= i__8; ++l) {
270 /*< ajb = a(jb+j) >*/
271 ajb = a[jb + j];
272 /*< aje = a(je+j) >*/
273 aje = a[je + j];
274 /*< t1 = ajb + aje >*/
275 t1 = ajb + aje;
276 /*< ajc = a(jc+j) >*/
277 ajc = a[jc + j];
278 /*< ajd = a(jd+j) >*/
279 ajd = a[jd + j];
280 /*< t2 = ajc + ajd >*/
281 t2 = ajc + ajd;
282 /*< t3 = ajb - aje >*/
283 t3 = ajb - aje;
284 /*< t4 = ajc - ajd >*/
285 t4 = ajc - ajd;
286 /*< t5 = t1 + t2 >*/
287 t5 = t1 + t2;
288 /*< t6 = c1 * ( t1 - t2 ) >*/
289 t6 = c1 * (t1 - t2);
290 /*< aja = a(ja+j) >*/
291 aja = a[ja + j];
292 /*< t7 = aja - 0.25 * t5 >*/
293 t7 = aja - t5 * (float).25;
294 /*< a(ja+j) = aja + t5 >*/
295 a[ja + j] = aja + t5;
296 /*< t8 = t7 + t6 >*/
297 t8 = t7 + t6;
298 /*< t9 = t7 - t6 >*/
299 t9 = t7 - t6;
300 /*< t10 = c3 * t3 - c2 * t4 >*/
301 t10 = c3 * t3 - c2 * t4;
302 /*< t11 = c2 * t3 + c3 * t4 >*/
303 t11 = c2 * t3 + c3 * t4;
304 /*< bjb = b(jb+j) >*/
305 bjb = b[jb + j];
306 /*< bje = b(je+j) >*/
307 bje = b[je + j];
308 /*< u1 = bjb + bje >*/
309 u1 = bjb + bje;
310 /*< bjc = b(jc+j) >*/
311 bjc = b[jc + j];
312 /*< bjd = b(jd+j) >*/
313 bjd = b[jd + j];
314 /*< u2 = bjc + bjd >*/
315 u2 = bjc + bjd;
316 /*< u3 = bjb - bje >*/
317 u3 = bjb - bje;
318 /*< u4 = bjc - bjd >*/
319 u4 = bjc - bjd;
320 /*< u5 = u1 + u2 >*/
321 u5 = u1 + u2;
322 /*< u6 = c1 * ( u1 - u2 ) >*/
323 u6 = c1 * (u1 - u2);
324 /*< bja = b(ja+j) >*/
325 bja = b[ja + j];
326 /*< u7 = bja - 0.25 * u5 >*/
327 u7 = bja - u5 * (float).25;
328 /*< b(ja+j) = bja + u5 >*/
329 b[ja + j] = bja + u5;
330 /*< u8 = u7 + u6 >*/
331 u8 = u7 + u6;
332 /*< u9 = u7 - u6 >*/
333 u9 = u7 - u6;
334 /*< u10 = c3 * u3 - c2 * u4 >*/
335 u10 = c3 * u3 - c2 * u4;
336 /*< u11 = c2 * u3 + c3 * u4 >*/
337 u11 = c2 * u3 + c3 * u4;
338 /*< a(jb+j) = t8 - u11 >*/
339 a[jb + j] = t8 - u11;
340 /*< b(jb+j) = u8 + t11 >*/
341 b[jb + j] = u8 + t11;
342 /*< a(je+j) = t8 + u11 >*/
343 a[je + j] = t8 + u11;
344 /*< b(je+j) = u8 - t11 >*/
345 b[je + j] = u8 - t11;
346 /*< a(jc+j) = t9 - u10 >*/
347 a[jc + j] = t9 - u10;
348 /*< b(jc+j) = u9 + t10 >*/
349 b[jc + j] = u9 + t10;
350 /*< a(jd+j) = t9 + u10 >*/
351 a[jd + j] = t9 + u10;
352 /*< b(jd+j) = u9 - t10 >*/
353 b[jd + j] = u9 - t10;
354 /*< j = j + jump >*/
355 j += *jump;
356 /*< 110 continue >*/
357 /* L110: */
358 }
359
360 /*< else >*/
361 } else {
362
363 /* dir$ ivdep,shortloop */
364 /*< do 130 l = 1 , nvex >*/
365 i__8 = nvex;
366 for (l = 1; l <= i__8; ++l) {
367 /*< ajb = a(jb+j) >*/
368 ajb = a[jb + j];
369 /*< aje = a(je+j) >*/
370 aje = a[je + j];
371 /*< t1 = ajb + aje >*/
372 t1 = ajb + aje;
373 /*< ajc = a(jc+j) >*/
374 ajc = a[jc + j];
375 /*< ajd = a(jd+j) >*/
376 ajd = a[jd + j];
377 /*< t2 = ajc + ajd >*/
378 t2 = ajc + ajd;
379 /*< t3 = ajb - aje >*/
380 t3 = ajb - aje;
381 /*< t4 = ajc - ajd >*/
382 t4 = ajc - ajd;
383 /*< t5 = t1 + t2 >*/
384 t5 = t1 + t2;
385 /*< t6 = c1 * ( t1 - t2 ) >*/
386 t6 = c1 * (t1 - t2);
387 /*< aja = a(ja+j) >*/
388 aja = a[ja + j];
389 /*< t7 = aja - 0.25 * t5 >*/
390 t7 = aja - t5 * (float).25;
391 /*< a(ja+j) = aja + t5 >*/
392 a[ja + j] = aja + t5;
393 /*< t8 = t7 + t6 >*/
394 t8 = t7 + t6;
395 /*< t9 = t7 - t6 >*/
396 t9 = t7 - t6;
397 /*< t10 = c3 * t3 - c2 * t4 >*/
398 t10 = c3 * t3 - c2 * t4;
399 /*< t11 = c2 * t3 + c3 * t4 >*/
400 t11 = c2 * t3 + c3 * t4;
401 /*< bjb = b(jb+j) >*/
402 bjb = b[jb + j];
403 /*< bje = b(je+j) >*/
404 bje = b[je + j];
405 /*< u1 = bjb + bje >*/
406 u1 = bjb + bje;
407 /*< bjc = b(jc+j) >*/
408 bjc = b[jc + j];
409 /*< bjd = b(jd+j) >*/
410 bjd = b[jd + j];
411 /*< u2 = bjc + bjd >*/
412 u2 = bjc + bjd;
413 /*< u3 = bjb - bje >*/
414 u3 = bjb - bje;
415 /*< u4 = bjc - bjd >*/
416 u4 = bjc - bjd;
417 /*< u5 = u1 + u2 >*/
418 u5 = u1 + u2;
419 /*< u6 = c1 * ( u1 - u2 ) >*/
420 u6 = c1 * (u1 - u2);
421 /*< bja = b(ja+j) >*/
422 bja = b[ja + j];
423 /*< u7 = bja - 0.25 * u5 >*/
424 u7 = bja - u5 * (float).25;
425 /*< b(ja+j) = bja + u5 >*/
426 b[ja + j] = bja + u5;
427 /*< u8 = u7 + u6 >*/
428 u8 = u7 + u6;
429 /*< u9 = u7 - u6 >*/
430 u9 = u7 - u6;
431 /*< u10 = c3 * u3 - c2 * u4 >*/
432 u10 = c3 * u3 - c2 * u4;
433 /*< u11 = c2 * u3 + c3 * u4 >*/
434 u11 = c2 * u3 + c3 * u4;
435 /*< a(jb+j) = co1*(t8-u11) - si1*(u8+t11) >*/
436 a[jb + j] = co1 * (t8 - u11) - si1 * (u8 +
437 t11);
438 /*< b(jb+j) = si1*(t8-u11) + co1*(u8+t11) >*/
439 b[jb + j] = si1 * (t8 - u11) + co1 * (u8 +
440 t11);
441 /*< a(je+j) = co4*(t8+u11) - si4*(u8-t11) >*/
442 a[je + j] = co4 * (t8 + u11) - si4 * (u8 -
443 t11);
444 /*< b(je+j) = si4*(t8+u11) + co4*(u8-t11) >*/
445 b[je + j] = si4 * (t8 + u11) + co4 * (u8 -
446 t11);
447 /*< a(jc+j) = co2*(t9-u10) - si2*(u9+t10) >*/
448 a[jc + j] = co2 * (t9 - u10) - si2 * (u9 +
449 t10);
450 /*< b(jc+j) = si2*(t9-u10) + co2*(u9+t10) >*/
451 b[jc + j] = si2 * (t9 - u10) + co2 * (u9 +
452 t10);
453 /*< a(jd+j) = co3*(t9+u10) - si3*(u9-t10) >*/
454 a[jd + j] = co3 * (t9 + u10) - si3 * (u9 -
455 t10);
456 /*< b(jd+j) = si3*(t9+u10) + co3*(u9-t10) >*/
457 b[jd + j] = si3 * (t9 + u10) + co3 * (u9 -
458 t10);
459 /*< j = j + jump >*/
460 j += *jump;
461 /*< 130 continue >*/
462 /* L130: */
463 }
464
465 /*< endif >*/
466 }
467
468 /* -----( end of loop across transforms ) */
469
470 /*< ja = ja + jstepx >*/
471 ja += jstepx;
472 /*< if (ja.lt.istart) ja = ja + ninc >*/
473 if (ja < istart) {
474 ja += ninc;
475 }
476 /*< 135 continue >*/
477 /* L135: */
478 }
479 /*< 140 continue >*/
480 /* L140: */
481 }
482 /* -----( end of loop along transforms ) */
483 /*< kk = kk + 2*la >*/
484 kk += la << 1;
485 /*< 150 continue >*/
486 /* L150: */
487 }
488 /* -----( end of loop on nonzero k ) */
489 /*< la = 5*la >*/
490 la *= 5;
491 /*< 160 continue >*/
492 /* L160: */
493 }
494 /* -----( end of loop on type I radix-5 passes) */
495
496 /*< if (n.eq.5) go to 490 >*/
497 if (*n == 5) {
498 goto L490;
499 }
500
501 /* loop on type II radix-5 passes */
502 /* ------------------------------ */
503 /*< 400 continue >*/
504 /* L400: */
505
506 /*< do 480 ipass = mh+1 , m >*/
507 i__2 = m;
508 for (ipass = mh + 1; ipass <= i__2; ++ipass) {
509 /*< jstep = (n*inc) / (5*la) >*/
510 jstep = *n * *inc / (la * 5);
511 /*< jstepl = jstep - ninc >*/
512 jstepl = jstep - ninc;
513 /*< laincl = la * ink - ninc >*/
514 laincl = la * ink - ninc;
515 /*< kk = 0 >*/
516 kk = 0;
517
518 /* loop on k */
519 /* --------- */
520 /*< do 470 k = 0 , jstep-ink , ink >*/
521 i__4 = jstep - ink;
522 i__3 = ink;
523 for (k = 0; i__3 < 0 ? k >= i__4 : k <= i__4; k += i__3) {
524
525 /*< if (k.gt.0) then >*/
526 if (k > 0) {
527 /*< co1 = trigs(kk+1) >*/
528 co1 = trigs[kk + 1];
529 /*< si1 = s*trigs(kk+2) >*/
530 si1 = s * trigs[kk + 2];
531 /*< co2 = trigs(2*kk+1) >*/
532 co2 = trigs[(kk << 1) + 1];
533 /*< si2 = s*trigs(2*kk+2) >*/
534 si2 = s * trigs[(kk << 1) + 2];
535 /*< co3 = trigs(3*kk+1) >*/
536 co3 = trigs[kk * 3 + 1];
537 /*< si3 = s*trigs(3*kk+2) >*/
538 si3 = s * trigs[kk * 3 + 2];
539 /*< co4 = trigs(4*kk+1) >*/
540 co4 = trigs[(kk << 2) + 1];
541 /*< si4 = s*trigs(4*kk+2) >*/
542 si4 = s * trigs[(kk << 2) + 2];
543 /*< endif >*/
544 }
545
546 /* double loop along first transform in block */
547 /* ------------------------------------------ */
548 /*< do 460 ll = k , (la-1)*ink , 5*jstep >*/
549 i__6 = (la - 1) * ink;
550 i__5 = jstep * 5;
551 for (ll = k; i__5 < 0 ? ll >= i__6 : ll <= i__6; ll += i__5) {
552
553 /*< do 450 jjj = ll , (n-1)*inc , 5*la*ink >*/
554 i__7 = (*n - 1) * *inc;
555 i__8 = la * 5 * ink;
556 for (jjj = ll; i__8 < 0 ? jjj >= i__7 : jjj <= i__7; jjj
557 += i__8) {
558 /*< ja = istart + jjj >*/
559 ja = istart + jjj;
560
561 /* "transverse" loop */
562 /* ----------------- */
563 /*< do 445 nu = 1 , inq >*/
564 i__9 = inq;
565 for (nu = 1; nu <= i__9; ++nu) {
566 /*< jb = ja + jstepl >*/
567 jb = ja + jstepl;
568 /*< if (jb.lt.istart) jb = jb + ninc >*/
569 if (jb < istart) {
570 jb += ninc;
571 }
572 /*< jc = jb + jstepl >*/
573 jc = jb + jstepl;
574 /*< if (jc.lt.istart) jc = jc + ninc >*/
575 if (jc < istart) {
576 jc += ninc;
577 }
578 /*< jd = jc + jstepl >*/
579 jd = jc + jstepl;
580 /*< if (jd.lt.istart) jd = jd + ninc >*/
581 if (jd < istart) {
582 jd += ninc;
583 }
584 /*< je = jd + jstepl >*/
585 je = jd + jstepl;
586 /*< if (je.lt.istart) je = je + ninc >*/
587 if (je < istart) {
588 je += ninc;
589 }
590 /*< jf = ja + laincl >*/
591 jf = ja + laincl;
592 /*< if (jf.lt.istart) jf = jf + ninc >*/
593 if (jf < istart) {
594 jf += ninc;
595 }
596 /*< jg = jf + jstepl >*/
597 jg = jf + jstepl;
598 /*< if (jg.lt.istart) jg = jg + ninc >*/
599 if (jg < istart) {
600 jg += ninc;
601 }
602 /*< jh = jg + jstepl >*/
603 jh = jg + jstepl;
604 /*< if (jh.lt.istart) jh = jh + ninc >*/
605 if (jh < istart) {
606 jh += ninc;
607 }
608 /*< ji = jh + jstepl >*/
609 ji = jh + jstepl;
610 /*< if (ji.lt.istart) ji = ji + ninc >*/
611 if (ji < istart) {
612 ji += ninc;
613 }
614 /*< jj = ji + jstepl >*/
615 jj = ji + jstepl;
616 /*< if (jj.lt.istart) jj = jj + ninc >*/
617 if (jj < istart) {
618 jj += ninc;
619 }
620 /*< jk = jf + laincl >*/
621 jk = jf + laincl;
622 /*< if (jk.lt.istart) jk = jk + ninc >*/
623 if (jk < istart) {
624 jk += ninc;
625 }
626 /*< jl = jk + jstepl >*/
627 jl = jk + jstepl;
628 /*< if (jl.lt.istart) jl = jl + ninc >*/
629 if (jl < istart) {
630 jl += ninc;
631 }
632 /*< jm = jl + jstepl >*/
633 jm = jl + jstepl;
634 /*< if (jm.lt.istart) jm = jm + ninc >*/
635 if (jm < istart) {
636 jm += ninc;
637 }
638 /*< jn = jm + jstepl >*/
639 jn = jm + jstepl;
640 /*< if (jn.lt.istart) jn = jn + ninc >*/
641 if (jn < istart) {
642 jn += ninc;
643 }
644 /*< jo = jn + jstepl >*/
645 jo = jn + jstepl;
646 /*< if (jo.lt.istart) jo = jo + ninc >*/
647 if (jo < istart) {
648 jo += ninc;
649 }
650 /*< jp = jk + laincl >*/
651 jp = jk + laincl;
652 /*< if (jp.lt.istart) jp = jp + ninc >*/
653 if (jp < istart) {
654 jp += ninc;
655 }
656 /*< jq = jp + jstepl >*/
657 jq = jp + jstepl;
658 /*< if (jq.lt.istart) jq = jq + ninc >*/
659 if (jq < istart) {
660 jq += ninc;
661 }
662 /*< jr = jq + jstepl >*/
663 jr = jq + jstepl;
664 /*< if (jr.lt.istart) jr = jr + ninc >*/
665 if (jr < istart) {
666 jr += ninc;
667 }
668 /*< js = jr + jstepl >*/
669 js = jr + jstepl;
670 /*< if (js.lt.istart) js = js + ninc >*/
671 if (js < istart) {
672 js += ninc;
673 }
674 /*< jt = js + jstepl >*/
675 jt = js + jstepl;
676 /*< if (jt.lt.istart) jt = jt + ninc >*/
677 if (jt < istart) {
678 jt += ninc;
679 }
680 /*< ju = jp + laincl >*/
681 ju = jp + laincl;
682 /*< if (ju.lt.istart) ju = ju + ninc >*/
683 if (ju < istart) {
684 ju += ninc;
685 }
686 /*< jv = ju + jstepl >*/
687 jv = ju + jstepl;
688 /*< if (jv.lt.istart) jv = jv + ninc >*/
689 if (jv < istart) {
690 jv += ninc;
691 }
692 /*< jw = jv + jstepl >*/
693 jw = jv + jstepl;
694 /*< if (jw.lt.istart) jw = jw + ninc >*/
695 if (jw < istart) {
696 jw += ninc;
697 }
698 /*< jx = jw + jstepl >*/
699 jx = jw + jstepl;
700 /*< if (jx.lt.istart) jx = jx + ninc >*/
701 if (jx < istart) {
702 jx += ninc;
703 }
704 /*< jy = jx + jstepl >*/
705 jy = jx + jstepl;
706 /*< if (jy.lt.istart) jy = jy + ninc >*/
707 if (jy < istart) {
708 jy += ninc;
709 }
710 /*< j = 0 >*/
711 j = 0;
712
713 /* loop across transforms */
714 /* ---------------------- */
715 /*< if (k.eq.0) then >*/
716 if (k == 0) {
717
718 /* dir$ ivdep, shortloop */
719 /*< do 410 l = 1 , nvex >*/
720 i__10 = nvex;
721 for (l = 1; l <= i__10; ++l) {
722 /*< ajb = a(jb+j) >*/
723 ajb = a[jb + j];
724 /*< aje = a(je+j) >*/
725 aje = a[je + j];
726 /*< t1 = ajb + aje >*/
727 t1 = ajb + aje;
728 /*< ajc = a(jc+j) >*/
729 ajc = a[jc + j];
730 /*< ajd = a(jd+j) >*/
731 ajd = a[jd + j];
732 /*< t2 = ajc + ajd >*/
733 t2 = ajc + ajd;
734 /*< t3 = ajb - aje >*/
735 t3 = ajb - aje;
736 /*< t4 = ajc - ajd >*/
737 t4 = ajc - ajd;
738 /*< ajf = a(jf+j) >*/
739 ajf = a[jf + j];
740 /*< ajb = ajf >*/
741 ajb = ajf;
742 /*< t5 = t1 + t2 >*/
743 t5 = t1 + t2;
744 /*< t6 = c1 * ( t1 - t2 ) >*/
745 t6 = c1 * (t1 - t2);
746 /*< aja = a(ja+j) >*/
747 aja = a[ja + j];
748 /*< t7 = aja - 0.25 * t5 >*/
749 t7 = aja - t5 * (float).25;
750 /*< a(ja+j) = aja + t5 >*/
751 a[ja + j] = aja + t5;
752 /*< t8 = t7 + t6 >*/
753 t8 = t7 + t6;
754 /*< t9 = t7 - t6 >*/
755 t9 = t7 - t6;
756 /*< ajk = a(jk+j) >*/
757 ajk = a[jk + j];
758 /*< ajc = ajk >*/
759 ajc = ajk;
760 /*< t10 = c3 * t3 - c2 * t4 >*/
761 t10 = c3 * t3 - c2 * t4;
762 /*< t11 = c2 * t3 + c3 * t4 >*/
763 t11 = c2 * t3 + c3 * t4;
764 /*< bjb = b(jb+j) >*/
765 bjb = b[jb + j];
766 /*< bje = b(je+j) >*/
767 bje = b[je + j];
768 /*< u1 = bjb + bje >*/
769 u1 = bjb + bje;
770 /*< bjc = b(jc+j) >*/
771 bjc = b[jc + j];
772 /*< bjd = b(jd+j) >*/
773 bjd = b[jd + j];
774 /*< u2 = bjc + bjd >*/
775 u2 = bjc + bjd;
776 /*< u3 = bjb - bje >*/
777 u3 = bjb - bje;
778 /*< u4 = bjc - bjd >*/
779 u4 = bjc - bjd;
780 /*< bjf = b(jf+j) >*/
781 bjf = b[jf + j];
782 /*< bjb = bjf >*/
783 bjb = bjf;
784 /*< u5 = u1 + u2 >*/
785 u5 = u1 + u2;
786 /*< u6 = c1 * ( u1 - u2 ) >*/
787 u6 = c1 * (u1 - u2);
788 /*< bja = b(ja+j) >*/
789 bja = b[ja + j];
790 /*< u7 = bja - 0.25 * u5 >*/
791 u7 = bja - u5 * (float).25;
792 /*< b(ja+j) = bja + u5 >*/
793 b[ja + j] = bja + u5;
794 /*< u8 = u7 + u6 >*/
795 u8 = u7 + u6;
796 /*< u9 = u7 - u6 >*/
797 u9 = u7 - u6;
798 /*< bjk = b(jk+j) >*/
799 bjk = b[jk + j];
800 /*< bjc = bjk >*/
801 bjc = bjk;
802 /*< u10 = c3 * u3 - c2 * u4 >*/
803 u10 = c3 * u3 - c2 * u4;
804 /*< u11 = c2 * u3 + c3 * u4 >*/
805 u11 = c2 * u3 + c3 * u4;
806 /*< a(jf+j) = t8 - u11 >*/
807 a[jf + j] = t8 - u11;
808 /*< b(jf+j) = u8 + t11 >*/
809 b[jf + j] = u8 + t11;
810 /*< aje = t8 + u11 >*/
811 aje = t8 + u11;
812 /*< bje = u8 - t11 >*/
813 bje = u8 - t11;
814 /*< a(jk+j) = t9 - u10 >*/
815 a[jk + j] = t9 - u10;
816 /*< b(jk+j) = u9 + t10 >*/
817 b[jk + j] = u9 + t10;
818 /*< ajd = t9 + u10 >*/
819 ajd = t9 + u10;
820 /*< bjd = u9 - t10 >*/
821 bjd = u9 - t10;
822 /* ---------------------- */
823 /*< ajg = a(jg+j) >*/
824 ajg = a[jg + j];
825 /*< ajj = a(jj+j) >*/
826 ajj = a[jj + j];
827 /*< t1 = ajg + ajj >*/
828 t1 = ajg + ajj;
829 /*< ajh = a(jh+j) >*/
830 ajh = a[jh + j];
831 /*< aji = a(ji+j) >*/
832 aji = a[ji + j];
833 /*< t2 = ajh + aji >*/
834 t2 = ajh + aji;
835 /*< t3 = ajg - ajj >*/
836 t3 = ajg - ajj;
837 /*< t4 = ajh - aji >*/
838 t4 = ajh - aji;
839 /*< ajl = a(jl+j) >*/
840 ajl = a[jl + j];
841 /*< ajh = ajl >*/
842 ajh = ajl;
843 /*< t5 = t1 + t2 >*/
844 t5 = t1 + t2;
845 /*< t6 = c1 * ( t1 - t2 ) >*/
846 t6 = c1 * (t1 - t2);
847 /*< t7 = ajb - 0.25 * t5 >*/
848 t7 = ajb - t5 * (float).25;
849 /*< a(jb+j) = ajb + t5 >*/
850 a[jb + j] = ajb + t5;
851 /*< t8 = t7 + t6 >*/
852 t8 = t7 + t6;
853 /*< t9 = t7 - t6 >*/
854 t9 = t7 - t6;
855 /*< ajq = a(jq+j) >*/
856 ajq = a[jq + j];
857 /*< aji = ajq >*/
858 aji = ajq;
859 /*< t10 = c3 * t3 - c2 * t4 >*/
860 t10 = c3 * t3 - c2 * t4;
861 /*< t11 = c2 * t3 + c3 * t4 >*/
862 t11 = c2 * t3 + c3 * t4;
863 /*< bjg = b(jg+j) >*/
864 bjg = b[jg + j];
865 /*< bjj = b(jj+j) >*/
866 bjj = b[jj + j];
867 /*< u1 = bjg + bjj >*/
868 u1 = bjg + bjj;
869 /*< bjh = b(jh+j) >*/
870 bjh = b[jh + j];
871 /*< bji = b(ji+j) >*/
872 bji = b[ji + j];
873 /*< u2 = bjh + bji >*/
874 u2 = bjh + bji;
875 /*< u3 = bjg - bjj >*/
876 u3 = bjg - bjj;
877 /*< u4 = bjh - bji >*/
878 u4 = bjh - bji;
879 /*< bjl = b(jl+j) >*/
880 bjl = b[jl + j];
881 /*< bjh = bjl >*/
882 bjh = bjl;
883 /*< u5 = u1 + u2 >*/
884 u5 = u1 + u2;
885 /*< u6 = c1 * ( u1 - u2 ) >*/
886 u6 = c1 * (u1 - u2);
887 /*< u7 = bjb - 0.25 * u5 >*/
888 u7 = bjb - u5 * (float).25;
889 /*< b(jb+j) = bjb + u5 >*/
890 b[jb + j] = bjb + u5;
891 /*< u8 = u7 + u6 >*/
892 u8 = u7 + u6;
893 /*< u9 = u7 - u6 >*/
894 u9 = u7 - u6;
895 /*< bjq = b(jq+j) >*/
896 bjq = b[jq + j];
897 /*< bji = bjq >*/
898 bji = bjq;
899 /*< u10 = c3 * u3 - c2 * u4 >*/
900 u10 = c3 * u3 - c2 * u4;
901 /*< u11 = c2 * u3 + c3 * u4 >*/
902 u11 = c2 * u3 + c3 * u4;
903 /*< a(jg+j) = t8 - u11 >*/
904 a[jg + j] = t8 - u11;
905 /*< b(jg+j) = u8 + t11 >*/
906 b[jg + j] = u8 + t11;
907 /*< ajj = t8 + u11 >*/
908 ajj = t8 + u11;
909 /*< bjj = u8 - t11 >*/
910 bjj = u8 - t11;
911 /*< a(jl+j) = t9 - u10 >*/
912 a[jl + j] = t9 - u10;
913 /*< b(jl+j) = u9 + t10 >*/
914 b[jl + j] = u9 + t10;
915 /*< a(jq+j) = t9 + u10 >*/
916 a[jq + j] = t9 + u10;
917 /*< b(jq+j) = u9 - t10 >*/
918 b[jq + j] = u9 - t10;
919 /* ---------------------- */
920 /*< ajo = a(jo+j) >*/
921 ajo = a[jo + j];
922 /*< t1 = ajh + ajo >*/
923 t1 = ajh + ajo;
924 /*< ajm = a(jm+j) >*/
925 ajm = a[jm + j];
926 /*< ajn = a(jn+j) >*/
927 ajn = a[jn + j];
928 /*< t2 = ajm + ajn >*/
929 t2 = ajm + ajn;
930 /*< t3 = ajh - ajo >*/
931 t3 = ajh - ajo;
932 /*< t4 = ajm - ajn >*/
933 t4 = ajm - ajn;
934 /*< ajr = a(jr+j) >*/
935 ajr = a[jr + j];
936 /*< ajn = ajr >*/
937 ajn = ajr;
938 /*< t5 = t1 + t2 >*/
939 t5 = t1 + t2;
940 /*< t6 = c1 * ( t1 - t2 ) >*/
941 t6 = c1 * (t1 - t2);
942 /*< t7 = ajc - 0.25 * t5 >*/
943 t7 = ajc - t5 * (float).25;
944 /*< a(jc+j) = ajc + t5 >*/
945 a[jc + j] = ajc + t5;
946 /*< t8 = t7 + t6 >*/
947 t8 = t7 + t6;
948 /*< t9 = t7 - t6 >*/
949 t9 = t7 - t6;
950 /*< ajw = a(jw+j) >*/
951 ajw = a[jw + j];
952 /*< ajo = ajw >*/
953 ajo = ajw;
954 /*< t10 = c3 * t3 - c2 * t4 >*/
955 t10 = c3 * t3 - c2 * t4;
956 /*< t11 = c2 * t3 + c3 * t4 >*/
957 t11 = c2 * t3 + c3 * t4;
958 /*< bjo = b(jo+j) >*/
959 bjo = b[jo + j];
960 /*< u1 = bjh + bjo >*/
961 u1 = bjh + bjo;
962 /*< bjm = b(jm+j) >*/
963 bjm = b[jm + j];
964 /*< bjn = b(jn+j) >*/
965 bjn = b[jn + j];
966 /*< u2 = bjm + bjn >*/
967 u2 = bjm + bjn;
968 /*< u3 = bjh - bjo >*/
969 u3 = bjh - bjo;
970 /*< u4 = bjm - bjn >*/
971 u4 = bjm - bjn;
972 /*< bjr = b(jr+j) >*/
973 bjr = b[jr + j];
974 /*< bjn = bjr >*/
975 bjn = bjr;
976 /*< u5 = u1 + u2 >*/
977 u5 = u1 + u2;
978 /*< u6 = c1 * ( u1 - u2 ) >*/
979 u6 = c1 * (u1 - u2);
980 /*< u7 = bjc - 0.25 * u5 >*/
981 u7 = bjc - u5 * (float).25;
982 /*< b(jc+j) = bjc + u5 >*/
983 b[jc + j] = bjc + u5;
984 /*< u8 = u7 + u6 >*/
985 u8 = u7 + u6;
986 /*< u9 = u7 - u6 >*/
987 u9 = u7 - u6;
988 /*< bjw = b(jw+j) >*/
989 bjw = b[jw + j];
990 /*< bjo = bjw >*/
991 bjo = bjw;
992 /*< u10 = c3 * u3 - c2 * u4 >*/
993 u10 = c3 * u3 - c2 * u4;
994 /*< u11 = c2 * u3 + c3 * u4 >*/
995 u11 = c2 * u3 + c3 * u4;
996 /*< a(jh+j) = t8 - u11 >*/
997 a[jh + j] = t8 - u11;
998 /*< b(jh+j) = u8 + t11 >*/
999 b[jh + j] = u8 + t11;
1000 /*< a(jw+j) = t8 + u11 >*/
1001 a[jw + j] = t8 + u11;
1002 /*< b(jw+j) = u8 - t11 >*/
1003 b[jw + j] = u8 - t11;
1004 /*< a(jm+j) = t9 - u10 >*/
1005 a[jm + j] = t9 - u10;
1006 /*< b(jm+j) = u9 + t10 >*/
1007 b[jm + j] = u9 + t10;
1008 /*< a(jr+j) = t9 + u10 >*/
1009 a[jr + j] = t9 + u10;
1010 /*< b(jr+j) = u9 - t10 >*/
1011 b[jr + j] = u9 - t10;
1012 /* ---------------------- */
1013 /*< ajt = a(jt+j) >*/
1014 ajt = a[jt + j];
1015 /*< t1 = aji + ajt >*/
1016 t1 = aji + ajt;
1017 /*< ajs = a(js+j) >*/
1018 ajs = a[js + j];
1019 /*< t2 = ajn + ajs >*/
1020 t2 = ajn + ajs;
1021 /*< t3 = aji - ajt >*/
1022 t3 = aji - ajt;
1023 /*< t4 = ajn - ajs >*/
1024 t4 = ajn - ajs;
1025 /*< ajx = a(jx+j) >*/
1026 ajx = a[jx + j];
1027 /*< ajt = ajx >*/
1028 ajt = ajx;
1029 /*< t5 = t1 + t2 >*/
1030 t5 = t1 + t2;
1031 /*< t6 = c1 * ( t1 - t2 ) >*/
1032 t6 = c1 * (t1 - t2);
1033 /*< ajp = a(jp+j) >*/
1034 ajp = a[jp + j];
1035 /*< t7 = ajp - 0.25 * t5 >*/
1036 t7 = ajp - t5 * (float).25;
1037 /*< ax = ajp + t5 >*/
1038 ax = ajp + t5;
1039 /*< t8 = t7 + t6 >*/
1040 t8 = t7 + t6;
1041 /*< t9 = t7 - t6 >*/
1042 t9 = t7 - t6;
1043 /*< a(jp+j) = ajd >*/
1044 a[jp + j] = ajd;
1045 /*< t10 = c3 * t3 - c2 * t4 >*/
1046 t10 = c3 * t3 - c2 * t4;
1047 /*< t11 = c2 * t3 + c3 * t4 >*/
1048 t11 = c2 * t3 + c3 * t4;
1049 /*< a(jd+j) = ax >*/
1050 a[jd + j] = ax;
1051 /*< bjt = b(jt+j) >*/
1052 bjt = b[jt + j];
1053 /*< u1 = bji + bjt >*/
1054 u1 = bji + bjt;
1055 /*< bjs = b(js+j) >*/
1056 bjs = b[js + j];
1057 /*< u2 = bjn + bjs >*/
1058 u2 = bjn + bjs;
1059 /*< u3 = bji - bjt >*/
1060 u3 = bji - bjt;
1061 /*< u4 = bjn - bjs >*/
1062 u4 = bjn - bjs;
1063 /*< bjx = b(jx+j) >*/
1064 bjx = b[jx + j];
1065 /*< bjt = bjx >*/
1066 bjt = bjx;
1067 /*< u5 = u1 + u2 >*/
1068 u5 = u1 + u2;
1069 /*< u6 = c1 * ( u1 - u2 ) >*/
1070 u6 = c1 * (u1 - u2);
1071 /*< bjp = b(jp+j) >*/
1072 bjp = b[jp + j];
1073 /*< u7 = bjp - 0.25 * u5 >*/
1074 u7 = bjp - u5 * (float).25;
1075 /*< bx = bjp + u5 >*/
1076 bx = bjp + u5;
1077 /*< u8 = u7 + u6 >*/
1078 u8 = u7 + u6;
1079 /*< u9 = u7 - u6 >*/
1080 u9 = u7 - u6;
1081 /*< b(jp+j) = bjd >*/
1082 b[jp + j] = bjd;
1083 /*< u10 = c3 * u3 - c2 * u4 >*/
1084 u10 = c3 * u3 - c2 * u4;
1085 /*< u11 = c2 * u3 + c3 * u4 >*/
1086 u11 = c2 * u3 + c3 * u4;
1087 /*< b(jd+j) = bx >*/
1088 b[jd + j] = bx;
1089 /*< a(ji+j) = t8 - u11 >*/
1090 a[ji + j] = t8 - u11;
1091 /*< b(ji+j) = u8 + t11 >*/
1092 b[ji + j] = u8 + t11;
1093 /*< a(jx+j) = t8 + u11 >*/
1094 a[jx + j] = t8 + u11;
1095 /*< b(jx+j) = u8 - t11 >*/
1096 b[jx + j] = u8 - t11;
1097 /*< a(jn+j) = t9 - u10 >*/
1098 a[jn + j] = t9 - u10;
1099 /*< b(jn+j) = u9 + t10 >*/
1100 b[jn + j] = u9 + t10;
1101 /*< a(js+j) = t9 + u10 >*/
1102 a[js + j] = t9 + u10;
1103 /*< b(js+j) = u9 - t10 >*/
1104 b[js + j] = u9 - t10;
1105 /* ---------------------- */
1106 /*< ajv = a(jv+j) >*/
1107 ajv = a[jv + j];
1108 /*< ajy = a(jy+j) >*/
1109 ajy = a[jy + j];
1110 /*< t1 = ajv + ajy >*/
1111 t1 = ajv + ajy;
1112 /*< t2 = ajo + ajt >*/
1113 t2 = ajo + ajt;
1114 /*< t3 = ajv - ajy >*/
1115 t3 = ajv - ajy;
1116 /*< t4 = ajo - ajt >*/
1117 t4 = ajo - ajt;
1118 /*< a(jv+j) = ajj >*/
1119 a[jv + j] = ajj;
1120 /*< t5 = t1 + t2 >*/
1121 t5 = t1 + t2;
1122 /*< t6 = c1 * ( t1 - t2 ) >*/
1123 t6 = c1 * (t1 - t2);
1124 /*< aju = a(ju+j) >*/
1125 aju = a[ju + j];
1126 /*< t7 = aju - 0.25 * t5 >*/
1127 t7 = aju - t5 * (float).25;
1128 /*< ax = aju + t5 >*/
1129 ax = aju + t5;
1130 /*< t8 = t7 + t6 >*/
1131 t8 = t7 + t6;
1132 /*< t9 = t7 - t6 >*/
1133 t9 = t7 - t6;
1134 /*< a(ju+j) = aje >*/
1135 a[ju + j] = aje;
1136 /*< t10 = c3 * t3 - c2 * t4 >*/
1137 t10 = c3 * t3 - c2 * t4;
1138 /*< t11 = c2 * t3 + c3 * t4 >*/
1139 t11 = c2 * t3 + c3 * t4;
1140 /*< a(je+j) = ax >*/
1141 a[je + j] = ax;
1142 /*< bjv = b(jv+j) >*/
1143 bjv = b[jv + j];
1144 /*< bjy = b(jy+j) >*/
1145 bjy = b[jy + j];
1146 /*< u1 = bjv + bjy >*/
1147 u1 = bjv + bjy;
1148 /*< u2 = bjo + bjt >*/
1149 u2 = bjo + bjt;
1150 /*< u3 = bjv - bjy >*/
1151 u3 = bjv - bjy;
1152 /*< u4 = bjo - bjt >*/
1153 u4 = bjo - bjt;
1154 /*< b(jv+j) = bjj >*/
1155 b[jv + j] = bjj;
1156 /*< u5 = u1 + u2 >*/
1157 u5 = u1 + u2;
1158 /*< u6 = c1 * ( u1 - u2 ) >*/
1159 u6 = c1 * (u1 - u2);
1160 /*< bju = b(ju+j) >*/
1161 bju = b[ju + j];
1162 /*< u7 = bju - 0.25 * u5 >*/
1163 u7 = bju - u5 * (float).25;
1164 /*< bx = bju + u5 >*/
1165 bx = bju + u5;
1166 /*< u8 = u7 + u6 >*/
1167 u8 = u7 + u6;
1168 /*< u9 = u7 - u6 >*/
1169 u9 = u7 - u6;
1170 /*< b(ju+j) = bje >*/
1171 b[ju + j] = bje;
1172 /*< u10 = c3 * u3 - c2 * u4 >*/
1173 u10 = c3 * u3 - c2 * u4;
1174 /*< u11 = c2 * u3 + c3 * u4 >*/
1175 u11 = c2 * u3 + c3 * u4;
1176 /*< b(je+j) = bx >*/
1177 b[je + j] = bx;
1178 /*< a(jj+j) = t8 - u11 >*/
1179 a[jj + j] = t8 - u11;
1180 /*< b(jj+j) = u8 + t11 >*/
1181 b[jj + j] = u8 + t11;
1182 /*< a(jy+j) = t8 + u11 >*/
1183 a[jy + j] = t8 + u11;
1184 /*< b(jy+j) = u8 - t11 >*/
1185 b[jy + j] = u8 - t11;
1186 /*< a(jo+j) = t9 - u10 >*/
1187 a[jo + j] = t9 - u10;
1188 /*< b(jo+j) = u9 + t10 >*/
1189 b[jo + j] = u9 + t10;
1190 /*< a(jt+j) = t9 + u10 >*/
1191 a[jt + j] = t9 + u10;
1192 /*< b(jt+j) = u9 - t10 >*/
1193 b[jt + j] = u9 - t10;
1194 /*< j = j + jump >*/
1195 j += *jump;
1196 /*< 410 continue >*/
1197 /* L410: */
1198 }
1199
1200 /*< else >*/
1201 } else {
1202
1203 /* dir$ ivdep, shortloop */
1204 /*< do 440 l = 1 , nvex >*/
1205 i__10 = nvex;
1206 for (l = 1; l <= i__10; ++l) {
1207 /*< ajb = a(jb+j) >*/
1208 ajb = a[jb + j];
1209 /*< aje = a(je+j) >*/
1210 aje = a[je + j];
1211 /*< t1 = ajb + aje >*/
1212 t1 = ajb + aje;
1213 /*< ajc = a(jc+j) >*/
1214 ajc = a[jc + j];
1215 /*< ajd = a(jd+j) >*/
1216 ajd = a[jd + j];
1217 /*< t2 = ajc + ajd >*/
1218 t2 = ajc + ajd;
1219 /*< t3 = ajb - aje >*/
1220 t3 = ajb - aje;
1221 /*< t4 = ajc - ajd >*/
1222 t4 = ajc - ajd;
1223 /*< ajf = a(jf+j) >*/
1224 ajf = a[jf + j];
1225 /*< ajb = ajf >*/
1226 ajb = ajf;
1227 /*< t5 = t1 + t2 >*/
1228 t5 = t1 + t2;
1229 /*< t6 = c1 * ( t1 - t2 ) >*/
1230 t6 = c1 * (t1 - t2);
1231 /*< aja = a(ja+j) >*/
1232 aja = a[ja + j];
1233 /*< t7 = aja - 0.25 * t5 >*/
1234 t7 = aja - t5 * (float).25;
1235 /*< a(ja+j) = aja + t5 >*/
1236 a[ja + j] = aja + t5;
1237 /*< t8 = t7 + t6 >*/
1238 t8 = t7 + t6;
1239 /*< t9 = t7 - t6 >*/
1240 t9 = t7 - t6;
1241 /*< ajk = a(jk+j) >*/
1242 ajk = a[jk + j];
1243 /*< ajc = ajk >*/
1244 ajc = ajk;
1245 /*< t10 = c3 * t3 - c2 * t4 >*/
1246 t10 = c3 * t3 - c2 * t4;
1247 /*< t11 = c2 * t3 + c3 * t4 >*/
1248 t11 = c2 * t3 + c3 * t4;
1249 /*< bjb = b(jb+j) >*/
1250 bjb = b[jb + j];
1251 /*< bje = b(je+j) >*/
1252 bje = b[je + j];
1253 /*< u1 = bjb + bje >*/
1254 u1 = bjb + bje;
1255 /*< bjc = b(jc+j) >*/
1256 bjc = b[jc + j];
1257 /*< bjd = b(jd+j) >*/
1258 bjd = b[jd + j];
1259 /*< u2 = bjc + bjd >*/
1260 u2 = bjc + bjd;
1261 /*< u3 = bjb - bje >*/
1262 u3 = bjb - bje;
1263 /*< u4 = bjc - bjd >*/
1264 u4 = bjc - bjd;
1265 /*< bjf = b(jf+j) >*/
1266 bjf = b[jf + j];
1267 /*< bjb = bjf >*/
1268 bjb = bjf;
1269 /*< u5 = u1 + u2 >*/
1270 u5 = u1 + u2;
1271 /*< u6 = c1 * ( u1 - u2 ) >*/
1272 u6 = c1 * (u1 - u2);
1273 /*< bja = b(ja+j) >*/
1274 bja = b[ja + j];
1275 /*< u7 = bja - 0.25 * u5 >*/
1276 u7 = bja - u5 * (float).25;
1277 /*< b(ja+j) = bja + u5 >*/
1278 b[ja + j] = bja + u5;
1279 /*< u8 = u7 + u6 >*/
1280 u8 = u7 + u6;
1281 /*< u9 = u7 - u6 >*/
1282 u9 = u7 - u6;
1283 /*< bjk = b(jk+j) >*/
1284 bjk = b[jk + j];
1285 /*< bjc = bjk >*/
1286 bjc = bjk;
1287 /*< u10 = c3 * u3 - c2 * u4 >*/
1288 u10 = c3 * u3 - c2 * u4;
1289 /*< u11 = c2 * u3 + c3 * u4 >*/
1290 u11 = c2 * u3 + c3 * u4;
1291 /*< a(jf+j) = co1*(t8-u11) - si1*(u8+t11) >*/
1292 a[jf + j] = co1 * (t8 - u11) - si1 * (u8
1293 + t11);
1294 /*< b(jf+j) = si1*(t8-u11) + co1*(u8+t11) >*/
1295 b[jf + j] = si1 * (t8 - u11) + co1 * (u8
1296 + t11);
1297 /*< aje = co4*(t8+u11) - si4*(u8-t11) >*/
1298 aje = co4 * (t8 + u11) - si4 * (u8 - t11);
1299 /*< bje = si4*(t8+u11) + co4*(u8-t11) >*/
1300 bje = si4 * (t8 + u11) + co4 * (u8 - t11);
1301 /*< a(jk+j) = co2*(t9-u10) - si2*(u9+t10) >*/
1302 a[jk + j] = co2 * (t9 - u10) - si2 * (u9
1303 + t10);
1304 /*< b(jk+j) = si2*(t9-u10) + co2*(u9+t10) >*/
1305 b[jk + j] = si2 * (t9 - u10) + co2 * (u9
1306 + t10);
1307 /*< ajd = co3*(t9+u10) - si3*(u9-t10) >*/
1308 ajd = co3 * (t9 + u10) - si3 * (u9 - t10);
1309 /*< bjd = si3*(t9+u10) + co3*(u9-t10) >*/
1310 bjd = si3 * (t9 + u10) + co3 * (u9 - t10);
1311 /* ---------------------- */
1312 /*< ajg = a(jg+j) >*/
1313 ajg = a[jg + j];
1314 /*< ajj = a(jj+j) >*/
1315 ajj = a[jj + j];
1316 /*< t1 = ajg + ajj >*/
1317 t1 = ajg + ajj;
1318 /*< ajh = a(jh+j) >*/
1319 ajh = a[jh + j];
1320 /*< aji = a(ji+j) >*/
1321 aji = a[ji + j];
1322 /*< t2 = ajh + aji >*/
1323 t2 = ajh + aji;
1324 /*< t3 = ajg - ajj >*/
1325 t3 = ajg - ajj;
1326 /*< t4 = ajh - aji >*/
1327 t4 = ajh - aji;
1328 /*< ajl = a(jl+j) >*/
1329 ajl = a[jl + j];
1330 /*< ajh = ajl >*/
1331 ajh = ajl;
1332 /*< t5 = t1 + t2 >*/
1333 t5 = t1 + t2;
1334 /*< t6 = c1 * ( t1 - t2 ) >*/
1335 t6 = c1 * (t1 - t2);
1336 /*< t7 = ajb - 0.25 * t5 >*/
1337 t7 = ajb - t5 * (float).25;
1338 /*< a(jb+j) = ajb + t5 >*/
1339 a[jb + j] = ajb + t5;
1340 /*< t8 = t7 + t6 >*/
1341 t8 = t7 + t6;
1342 /*< t9 = t7 - t6 >*/
1343 t9 = t7 - t6;
1344 /*< ajq = a(jq+j) >*/
1345 ajq = a[jq + j];
1346 /*< aji = ajq >*/
1347 aji = ajq;
1348 /*< t10 = c3 * t3 - c2 * t4 >*/
1349 t10 = c3 * t3 - c2 * t4;
1350 /*< t11 = c2 * t3 + c3 * t4 >*/
1351 t11 = c2 * t3 + c3 * t4;
1352 /*< bjg = b(jg+j) >*/
1353 bjg = b[jg + j];
1354 /*< bjj = b(jj+j) >*/
1355 bjj = b[jj + j];
1356 /*< u1 = bjg + bjj >*/
1357 u1 = bjg + bjj;
1358 /*< bjh = b(jh+j) >*/
1359 bjh = b[jh + j];
1360 /*< bji = b(ji+j) >*/
1361 bji = b[ji + j];
1362 /*< u2 = bjh + bji >*/
1363 u2 = bjh + bji;
1364 /*< u3 = bjg - bjj >*/
1365 u3 = bjg - bjj;
1366 /*< u4 = bjh - bji >*/
1367 u4 = bjh - bji;
1368 /*< bjl = b(jl+j) >*/
1369 bjl = b[jl + j];
1370 /*< bjh = bjl >*/
1371 bjh = bjl;
1372 /*< u5 = u1 + u2 >*/
1373 u5 = u1 + u2;
1374 /*< u6 = c1 * ( u1 - u2 ) >*/
1375 u6 = c1 * (u1 - u2);
1376 /*< u7 = bjb - 0.25 * u5 >*/
1377 u7 = bjb - u5 * (float).25;
1378 /*< b(jb+j) = bjb + u5 >*/
1379 b[jb + j] = bjb + u5;
1380 /*< u8 = u7 + u6 >*/
1381 u8 = u7 + u6;
1382 /*< u9 = u7 - u6 >*/
1383 u9 = u7 - u6;
1384 /*< bjq = b(jq+j) >*/
1385 bjq = b[jq + j];
1386 /*< bji = bjq >*/
1387 bji = bjq;
1388 /*< u10 = c3 * u3 - c2 * u4 >*/
1389 u10 = c3 * u3 - c2 * u4;
1390 /*< u11 = c2 * u3 + c3 * u4 >*/
1391 u11 = c2 * u3 + c3 * u4;
1392 /*< a(jg+j) = co1*(t8-u11) - si1*(u8+t11) >*/
1393 a[jg + j] = co1 * (t8 - u11) - si1 * (u8
1394 + t11);
1395 /*< b(jg+j) = si1*(t8-u11) + co1*(u8+t11) >*/
1396 b[jg + j] = si1 * (t8 - u11) + co1 * (u8
1397 + t11);
1398 /*< ajj = co4*(t8+u11) - si4*(u8-t11) >*/
1399 ajj = co4 * (t8 + u11) - si4 * (u8 - t11);
1400 /*< bjj = si4*(t8+u11) + co4*(u8-t11) >*/
1401 bjj = si4 * (t8 + u11) + co4 * (u8 - t11);
1402 /*< a(jl+j) = co2*(t9-u10) - si2*(u9+t10) >*/
1403 a[jl + j] = co2 * (t9 - u10) - si2 * (u9
1404 + t10);
1405 /*< b(jl+j) = si2*(t9-u10) + co2*(u9+t10) >*/
1406 b[jl + j] = si2 * (t9 - u10) + co2 * (u9
1407 + t10);
1408 /*< a(jq+j) = co3*(t9+u10) - si3*(u9-t10) >*/
1409 a[jq + j] = co3 * (t9 + u10) - si3 * (u9
1410 - t10);
1411 /*< b(jq+j) = si3*(t9+u10) + co3*(u9-t10) >*/
1412 b[jq + j] = si3 * (t9 + u10) + co3 * (u9
1413 - t10);
1414 /* ---------------------- */
1415 /*< ajo = a(jo+j) >*/
1416 ajo = a[jo + j];
1417 /*< t1 = ajh + ajo >*/
1418 t1 = ajh + ajo;
1419 /*< ajm = a(jm+j) >*/
1420 ajm = a[jm + j];
1421 /*< ajn = a(jn+j) >*/
1422 ajn = a[jn + j];
1423 /*< t2 = ajm + ajn >*/
1424 t2 = ajm + ajn;
1425 /*< t3 = ajh - ajo >*/
1426 t3 = ajh - ajo;
1427 /*< t4 = ajm - ajn >*/
1428 t4 = ajm - ajn;
1429 /*< ajr = a(jr+j) >*/
1430 ajr = a[jr + j];
1431 /*< ajn = ajr >*/
1432 ajn = ajr;
1433 /*< t5 = t1 + t2 >*/
1434 t5 = t1 + t2;
1435 /*< t6 = c1 * ( t1 - t2 ) >*/
1436 t6 = c1 * (t1 - t2);
1437 /*< t7 = ajc - 0.25 * t5 >*/
1438 t7 = ajc - t5 * (float).25;
1439 /*< a(jc+j) = ajc + t5 >*/
1440 a[jc + j] = ajc + t5;
1441 /*< t8 = t7 + t6 >*/
1442 t8 = t7 + t6;
1443 /*< t9 = t7 - t6 >*/
1444 t9 = t7 - t6;
1445 /*< ajw = a(jw+j) >*/
1446 ajw = a[jw + j];
1447 /*< ajo = ajw >*/
1448 ajo = ajw;
1449 /*< t10 = c3 * t3 - c2 * t4 >*/
1450 t10 = c3 * t3 - c2 * t4;
1451 /*< t11 = c2 * t3 + c3 * t4 >*/
1452 t11 = c2 * t3 + c3 * t4;
1453 /*< bjo = b(jo+j) >*/
1454 bjo = b[jo + j];
1455 /*< u1 = bjh + bjo >*/
1456 u1 = bjh + bjo;
1457 /*< bjm = b(jm+j) >*/
1458 bjm = b[jm + j];
1459 /*< bjn = b(jn+j) >*/
1460 bjn = b[jn + j];
1461 /*< u2 = bjm + bjn >*/
1462 u2 = bjm + bjn;
1463 /*< u3 = bjh - bjo >*/
1464 u3 = bjh - bjo;
1465 /*< u4 = bjm - bjn >*/
1466 u4 = bjm - bjn;
1467 /*< bjr = b(jr+j) >*/
1468 bjr = b[jr + j];
1469 /*< bjn = bjr >*/
1470 bjn = bjr;
1471 /*< u5 = u1 + u2 >*/
1472 u5 = u1 + u2;
1473 /*< u6 = c1 * ( u1 - u2 ) >*/
1474 u6 = c1 * (u1 - u2);
1475 /*< u7 = bjc - 0.25 * u5 >*/
1476 u7 = bjc - u5 * (float).25;
1477 /*< b(jc+j) = bjc + u5 >*/
1478 b[jc + j] = bjc + u5;
1479 /*< u8 = u7 + u6 >*/
1480 u8 = u7 + u6;
1481 /*< u9 = u7 - u6 >*/
1482 u9 = u7 - u6;
1483 /*< bjw = b(jw+j) >*/
1484 bjw = b[jw + j];
1485 /*< bjo = bjw >*/
1486 bjo = bjw;
1487 /*< u10 = c3 * u3 - c2 * u4 >*/
1488 u10 = c3 * u3 - c2 * u4;
1489 /*< u11 = c2 * u3 + c3 * u4 >*/
1490 u11 = c2 * u3 + c3 * u4;
1491 /*< a(jh+j) = co1*(t8-u11) - si1*(u8+t11) >*/
1492 a[jh + j] = co1 * (t8 - u11) - si1 * (u8
1493 + t11);
1494 /*< b(jh+j) = si1*(t8-u11) + co1*(u8+t11) >*/
1495 b[jh + j] = si1 * (t8 - u11) + co1 * (u8
1496 + t11);
1497 /*< a(jw+j) = co4*(t8+u11) - si4*(u8-t11) >*/
1498 a[jw + j] = co4 * (t8 + u11) - si4 * (u8
1499 - t11);
1500 /*< b(jw+j) = si4*(t8+u11) + co4*(u8-t11) >*/
1501 b[jw + j] = si4 * (t8 + u11) + co4 * (u8
1502 - t11);
1503 /*< a(jm+j) = co2*(t9-u10) - si2*(u9+t10) >*/
1504 a[jm + j] = co2 * (t9 - u10) - si2 * (u9
1505 + t10);
1506 /*< b(jm+j) = si2*(t9-u10) + co2*(u9+t10) >*/
1507 b[jm + j] = si2 * (t9 - u10) + co2 * (u9
1508 + t10);
1509 /*< a(jr+j) = co3*(t9+u10) - si3*(u9-t10) >*/
1510 a[jr + j] = co3 * (t9 + u10) - si3 * (u9
1511 - t10);
1512 /*< b(jr+j) = si3*(t9+u10) + co3*(u9-t10) >*/
1513 b[jr + j] = si3 * (t9 + u10) + co3 * (u9
1514 - t10);
1515 /* ---------------------- */
1516 /*< ajt = a(jt+j) >*/
1517 ajt = a[jt + j];
1518 /*< t1 = aji + ajt >*/
1519 t1 = aji + ajt;
1520 /*< ajs = a(js+j) >*/
1521 ajs = a[js + j];
1522 /*< t2 = ajn + ajs >*/
1523 t2 = ajn + ajs;
1524 /*< t3 = aji - ajt >*/
1525 t3 = aji - ajt;
1526 /*< t4 = ajn - ajs >*/
1527 t4 = ajn - ajs;
1528 /*< ajx = a(jx+j) >*/
1529 ajx = a[jx + j];
1530 /*< ajt = ajx >*/
1531 ajt = ajx;
1532 /*< t5 = t1 + t2 >*/
1533 t5 = t1 + t2;
1534 /*< t6 = c1 * ( t1 - t2 ) >*/
1535 t6 = c1 * (t1 - t2);
1536 /*< ajp = a(jp+j) >*/
1537 ajp = a[jp + j];
1538 /*< t7 = ajp - 0.25 * t5 >*/
1539 t7 = ajp - t5 * (float).25;
1540 /*< ax = ajp + t5 >*/
1541 ax = ajp + t5;
1542 /*< t8 = t7 + t6 >*/
1543 t8 = t7 + t6;
1544 /*< t9 = t7 - t6 >*/
1545 t9 = t7 - t6;
1546 /*< a(jp+j) = ajd >*/
1547 a[jp + j] = ajd;
1548 /*< t10 = c3 * t3 - c2 * t4 >*/
1549 t10 = c3 * t3 - c2 * t4;
1550 /*< t11 = c2 * t3 + c3 * t4 >*/
1551 t11 = c2 * t3 + c3 * t4;
1552 /*< a(jd+j) = ax >*/
1553 a[jd + j] = ax;
1554 /*< bjt = b(jt+j) >*/
1555 bjt = b[jt + j];
1556 /*< u1 = bji + bjt >*/
1557 u1 = bji + bjt;
1558 /*< bjs = b(js+j) >*/
1559 bjs = b[js + j];
1560 /*< u2 = bjn + bjs >*/
1561 u2 = bjn + bjs;
1562 /*< u3 = bji - bjt >*/
1563 u3 = bji - bjt;
1564 /*< u4 = bjn - bjs >*/
1565 u4 = bjn - bjs;
1566 /*< bjx = b(jx+j) >*/
1567 bjx = b[jx + j];
1568 /*< bjt = bjx >*/
1569 bjt = bjx;
1570 /*< u5 = u1 + u2 >*/
1571 u5 = u1 + u2;
1572 /*< u6 = c1 * ( u1 - u2 ) >*/
1573 u6 = c1 * (u1 - u2);
1574 /*< bjp = b(jp+j) >*/
1575 bjp = b[jp + j];
1576 /*< u7 = bjp - 0.25 * u5 >*/
1577 u7 = bjp - u5 * (float).25;
1578 /*< bx = bjp + u5 >*/
1579 bx = bjp + u5;
1580 /*< u8 = u7 + u6 >*/
1581 u8 = u7 + u6;
1582 /*< u9 = u7 - u6 >*/
1583 u9 = u7 - u6;
1584 /*< b(jp+j) = bjd >*/
1585 b[jp + j] = bjd;
1586 /*< u10 = c3 * u3 - c2 * u4 >*/
1587 u10 = c3 * u3 - c2 * u4;
1588 /*< u11 = c2 * u3 + c3 * u4 >*/
1589 u11 = c2 * u3 + c3 * u4;
1590 /*< b(jd+j) = bx >*/
1591 b[jd + j] = bx;
1592 /*< a(ji+j) = co1*(t8-u11) - si1*(u8+t11) >*/
1593 a[ji + j] = co1 * (t8 - u11) - si1 * (u8
1594 + t11);
1595 /*< b(ji+j) = si1*(t8-u11) + co1*(u8+t11) >*/
1596 b[ji + j] = si1 * (t8 - u11) + co1 * (u8
1597 + t11);
1598 /*< a(jx+j) = co4*(t8+u11) - si4*(u8-t11) >*/
1599 a[jx + j] = co4 * (t8 + u11) - si4 * (u8
1600 - t11);
1601 /*< b(jx+j) = si4*(t8+u11) + co4*(u8-t11) >*/
1602 b[jx + j] = si4 * (t8 + u11) + co4 * (u8
1603 - t11);
1604 /*< a(jn+j) = co2*(t9-u10) - si2*(u9+t10) >*/
1605 a[jn + j] = co2 * (t9 - u10) - si2 * (u9
1606 + t10);
1607 /*< b(jn+j) = si2*(t9-u10) + co2*(u9+t10) >*/
1608 b[jn + j] = si2 * (t9 - u10) + co2 * (u9
1609 + t10);
1610 /*< a(js+j) = co3*(t9+u10) - si3*(u9-t10) >*/
1611 a[js + j] = co3 * (t9 + u10) - si3 * (u9
1612 - t10);
1613 /*< b(js+j) = si3*(t9+u10) + co3*(u9-t10) >*/
1614 b[js + j] = si3 * (t9 + u10) + co3 * (u9
1615 - t10);
1616 /* ---------------------- */
1617 /*< ajv = a(jv+j) >*/
1618 ajv = a[jv + j];
1619 /*< ajy = a(jy+j) >*/
1620 ajy = a[jy + j];
1621 /*< t1 = ajv + ajy >*/
1622 t1 = ajv + ajy;
1623 /*< t2 = ajo + ajt >*/
1624 t2 = ajo + ajt;
1625 /*< t3 = ajv - ajy >*/
1626 t3 = ajv - ajy;
1627 /*< t4 = ajo - ajt >*/
1628 t4 = ajo - ajt;
1629 /*< a(jv+j) = ajj >*/
1630 a[jv + j] = ajj;
1631 /*< t5 = t1 + t2 >*/
1632 t5 = t1 + t2;
1633 /*< t6 = c1 * ( t1 - t2 ) >*/
1634 t6 = c1 * (t1 - t2);
1635 /*< aju = a(ju+j) >*/
1636 aju = a[ju + j];
1637 /*< t7 = aju - 0.25 * t5 >*/
1638 t7 = aju - t5 * (float).25;
1639 /*< ax = aju + t5 >*/
1640 ax = aju + t5;
1641 /*< t8 = t7 + t6 >*/
1642 t8 = t7 + t6;
1643 /*< t9 = t7 - t6 >*/
1644 t9 = t7 - t6;
1645 /*< a(ju+j) = aje >*/
1646 a[ju + j] = aje;
1647 /*< t10 = c3 * t3 - c2 * t4 >*/
1648 t10 = c3 * t3 - c2 * t4;
1649 /*< t11 = c2 * t3 + c3 * t4 >*/
1650 t11 = c2 * t3 + c3 * t4;
1651 /*< a(je+j) = ax >*/
1652 a[je + j] = ax;
1653 /*< bjv = b(jv+j) >*/
1654 bjv = b[jv + j];
1655 /*< bjy = b(jy+j) >*/
1656 bjy = b[jy + j];
1657 /*< u1 = bjv + bjy >*/
1658 u1 = bjv + bjy;
1659 /*< u2 = bjo + bjt >*/
1660 u2 = bjo + bjt;
1661 /*< u3 = bjv - bjy >*/
1662 u3 = bjv - bjy;
1663 /*< u4 = bjo - bjt >*/
1664 u4 = bjo - bjt;
1665 /*< b(jv+j) = bjj >*/
1666 b[jv + j] = bjj;
1667 /*< u5 = u1 + u2 >*/
1668 u5 = u1 + u2;
1669 /*< u6 = c1 * ( u1 - u2 ) >*/
1670 u6 = c1 * (u1 - u2);
1671 /*< bju = b(ju+j) >*/
1672 bju = b[ju + j];
1673 /*< u7 = bju - 0.25 * u5 >*/
1674 u7 = bju - u5 * (float).25;
1675 /*< bx = bju + u5 >*/
1676 bx = bju + u5;
1677 /*< u8 = u7 + u6 >*/
1678 u8 = u7 + u6;
1679 /*< u9 = u7 - u6 >*/
1680 u9 = u7 - u6;
1681 /*< b(ju+j) = bje >*/
1682 b[ju + j] = bje;
1683 /*< u10 = c3 * u3 - c2 * u4 >*/
1684 u10 = c3 * u3 - c2 * u4;
1685 /*< u11 = c2 * u3 + c3 * u4 >*/
1686 u11 = c2 * u3 + c3 * u4;
1687 /*< b(je+j) = bx >*/
1688 b[je + j] = bx;
1689 /*< a(jj+j) = co1*(t8-u11) - si1*(u8+t11) >*/
1690 a[jj + j] = co1 * (t8 - u11) - si1 * (u8
1691 + t11);
1692 /*< b(jj+j) = si1*(t8-u11) + co1*(u8+t11) >*/
1693 b[jj + j] = si1 * (t8 - u11) + co1 * (u8
1694 + t11);
1695 /*< a(jy+j) = co4*(t8+u11) - si4*(u8-t11) >*/
1696 a[jy + j] = co4 * (t8 + u11) - si4 * (u8
1697 - t11);
1698 /*< b(jy+j) = si4*(t8+u11) + co4*(u8-t11) >*/
1699 b[jy + j] = si4 * (t8 + u11) + co4 * (u8
1700 - t11);
1701 /*< a(jo+j) = co2*(t9-u10) - si2*(u9+t10) >*/
1702 a[jo + j] = co2 * (t9 - u10) - si2 * (u9
1703 + t10);
1704 /*< b(jo+j) = si2*(t9-u10) + co2*(u9+t10) >*/
1705 b[jo + j] = si2 * (t9 - u10) + co2 * (u9
1706 + t10);
1707 /*< a(jt+j) = co3*(t9+u10) - si3*(u9-t10) >*/
1708 a[jt + j] = co3 * (t9 + u10) - si3 * (u9
1709 - t10);
1710 /*< b(jt+j) = si3*(t9+u10) + co3*(u9-t10) >*/
1711 b[jt + j] = si3 * (t9 + u10) + co3 * (u9
1712 - t10);
1713 /*< j = j + jump >*/
1714 j += *jump;
1715 /*< 440 continue >*/
1716 /* L440: */
1717 }
1718
1719 /*< endif >*/
1720 }
1721
1722 /* -----(end of loop across transforms) */
1723
1724 /*< ja = ja + jstepx >*/
1725 ja += jstepx;
1726 /*< if (ja.lt.istart) ja = ja + ninc >*/
1727 if (ja < istart) {
1728 ja += ninc;
1729 }
1730 /*< 445 continue >*/
1731 /* L445: */
1732 }
1733 /*< 450 continue >*/
1734 /* L450: */
1735 }
1736 /*< 460 continue >*/
1737 /* L460: */
1738 }
1739 /* -----( end of double loop for this k ) */
1740 /*< kk = kk + 2*la >*/
1741 kk += la << 1;
1742 /*< 470 continue >*/
1743 /* L470: */
1744 }
1745 /* -----( end of loop over values of k ) */
1746 /*< la = 5*la >*/
1747 la *= 5;
1748 /*< 480 continue >*/
1749 /* L480: */
1750 }
1751 /* -----( end of loop on type II radix-5 passes ) */
1752 /* -----( nvex transforms completed) */
1753 /*< 490 continue >*/
1754 L490:
1755 /*< istart = istart + nvex * jump >*/
1756 istart += nvex * *jump;
1757 /*< 500 continue >*/
1758 /* L500: */
1759 }
1760 /* -----( end of loop on blocks of transforms ) */
1761
1762 /*< return >*/
1763 return 0;
1764 /*< end >*/
1765 } /* gpfa5f_ */
1766
1767 #ifdef __cplusplus
1768 }
1769 #endif
1770