Mercurial > octave-nkf

annotate libcruft/fftpack/cfftf1.f @ 7789:82be108cc558

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
First attempt at single precision tyeps * * * corrections to qrupdate single precision routines * * * prefer demotion to single over promotion to double * * * Add single precision support to log2 function * * * Trivial PROJECT file update * * * Cache optimized hermitian/transpose methods * * * Add tests for tranpose/hermitian and ChangeLog entry for new transpose code
author David Bateman <dbateman@free.fr>
date Sun, 27 Apr 2008 22:34:17 +0200
parents 44ed237bdc1e
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
1644
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
1 subroutine cfftf1 (n,c,ch,wa,ifac)
1645
44ed237bdc1e [project @ 1995-12-14 08:32:49 by jwe]
jwe
parents: 1644
diff changeset
2 dimension ch(*) ,c(*) ,wa(*) ,ifac(*)
1644
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
3 nf = ifac(2)
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
4 na = 0
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
5 l1 = 1
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
6 iw = 1
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
7 do 116 k1=1,nf
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
8 ip = ifac(k1+2)
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
9 l2 = ip*l1
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
10 ido = n/l2
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
11 idot = ido+ido
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
12 idl1 = idot*l1
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
13 if (ip .ne. 4) go to 103
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
14 ix2 = iw+idot
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
15 ix3 = ix2+idot
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
16 if (na .ne. 0) go to 101
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
17 call passf4 (idot,l1,c,ch,wa(iw),wa(ix2),wa(ix3))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
18 go to 102
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
19 101 call passf4 (idot,l1,ch,c,wa(iw),wa(ix2),wa(ix3))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
20 102 na = 1-na
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
21 go to 115
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
22 103 if (ip .ne. 2) go to 106
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
23 if (na .ne. 0) go to 104
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
24 call passf2 (idot,l1,c,ch,wa(iw))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
25 go to 105
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
26 104 call passf2 (idot,l1,ch,c,wa(iw))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
27 105 na = 1-na
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
28 go to 115
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
29 106 if (ip .ne. 3) go to 109
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
30 ix2 = iw+idot
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
31 if (na .ne. 0) go to 107
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
32 call passf3 (idot,l1,c,ch,wa(iw),wa(ix2))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
33 go to 108
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
34 107 call passf3 (idot,l1,ch,c,wa(iw),wa(ix2))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
35 108 na = 1-na
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
36 go to 115
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
37 109 if (ip .ne. 5) go to 112
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
38 ix2 = iw+idot
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
39 ix3 = ix2+idot
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
40 ix4 = ix3+idot
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
41 if (na .ne. 0) go to 110
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
42 call passf5 (idot,l1,c,ch,wa(iw),wa(ix2),wa(ix3),wa(ix4))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
43 go to 111
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
44 110 call passf5 (idot,l1,ch,c,wa(iw),wa(ix2),wa(ix3),wa(ix4))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
45 111 na = 1-na
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
46 go to 115
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
47 112 if (na .ne. 0) go to 113
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
48 call passf (nac,idot,ip,l1,idl1,c,c,c,ch,ch,wa(iw))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
49 go to 114
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
50 113 call passf (nac,idot,ip,l1,idl1,ch,ch,ch,c,c,wa(iw))
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
51 114 if (nac .ne. 0) na = 1-na
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
52 115 l1 = l2
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
53 iw = iw+(ip-1)*idot
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
54 116 continue
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
55 if (na .eq. 0) return
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
56 n2 = n+n
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
57 do 117 i=1,n2
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
58 c(i) = ch(i)
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
59 117 continue
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
60 return
395bb6d6c096 [project @ 1995-12-14 08:32:12 by jwe]
jwe
parents:
diff changeset
61 end