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