diff --git a/fourier/cfft.go b/fourier/cfft.go index f383a58d..0ba2d841 100644 --- a/fourier/cfft.go +++ b/fourier/cfft.go @@ -443,13 +443,12 @@ func passf(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) ch2m := newTwoArray(idl1, ip, ch2) idot := ido / 2 - ipp2 := ip + 1 ipph := (ip + 1) / 2 idp := ip * ido if ido < l1 { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for i := 0; i < ido; i++ { for k := 0; k < l1; k++ { ch3.set(i, k, j, cc3.at(i, j, k)+cc3.at(i, jc, k)) @@ -464,7 +463,7 @@ func passf(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) } } else { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for k := 0; k < l1; k++ { for i := 0; i < ido; i++ { ch3.set(i, k, j, cc3.at(i, j, k)+cc3.at(i, jc, k)) @@ -482,7 +481,7 @@ func passf(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) idl := 1 - ido inc := 0 for l := 1; l < ipph; l++ { - lc := ipp2 - (l + 1) + lc := ip - l idl += ido for ik := 0; ik < idl1; ik++ { c2m.set(ik, l, ch2m.at(ik, 0)+wa[idl-1]*ch2m.at(ik, 1)) @@ -491,7 +490,7 @@ func passf(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) idlj := idl inc += ido for j := 2; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j idlj += inc if idlj > idp { idlj -= idp @@ -512,7 +511,7 @@ func passf(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) } for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for ik := 1; ik < idl1; ik += 2 { ch2m.set(ik-1, j, c2m.at(ik-1, j)-c2m.at(ik, jc)) ch2m.set(ik-1, jc, c2m.at(ik-1, j)+c2m.at(ik, jc)) @@ -907,13 +906,12 @@ func passb(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) ch2m := newTwoArray(idl1, ip, ch2) idot := ido / 2 - ipp2 := ip + 1 ipph := (ip + 1) / 2 idp := ip * ido if ido < l1 { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for i := 0; i < ido; i++ { for k := 0; k < l1; k++ { ch3.set(i, k, j, cc3.at(i, j, k)+cc3.at(i, jc, k)) @@ -928,7 +926,7 @@ func passb(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) } } else { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for k := 0; k < l1; k++ { for i := 0; i < ido; i++ { ch3.set(i, k, j, cc3.at(i, j, k)+cc3.at(i, jc, k)) @@ -946,7 +944,7 @@ func passb(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) idl := 1 - ido inc := 0 for l := 1; l < ipph; l++ { - lc := ipp2 - (l + 1) + lc := ip - l idl += ido for ik := 0; ik < idl1; ik++ { c2m.set(ik, l, ch2m.at(ik, 0)+wa[idl-1]*ch2m.at(ik, 1)) @@ -955,7 +953,7 @@ func passb(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) idlj := idl inc += ido for j := 2; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j idlj += inc if idlj > idp { idlj -= idp @@ -976,7 +974,7 @@ func passb(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) (nac bool) } for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for ik := 1; ik < idl1; ik += 2 { ch2m.set(ik-1, j, c2m.at(ik-1, j)-c2m.at(ik, jc)) ch2m.set(ik-1, jc, c2m.at(ik-1, j)+c2m.at(ik, jc)) diff --git a/fourier/rfft.go b/fourier/rfft.go index 32d2d2ab..b3acc104 100644 --- a/fourier/rfft.go +++ b/fourier/rfft.go @@ -445,8 +445,6 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { dcp := math.Cos(arg) dsp := math.Sin(arg) ipph := (ip + 1) / 2 - ipp2 := ip + 1 - idp2 := ido + 1 nbd := (ido - 1) / 2 if ido == 1 { @@ -491,7 +489,7 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } if nbd < l1 { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for i := 2; i < ido; i += 2 { for k := 0; k < l1; k++ { c13.set(i-1, k, j, ch3.at(i-1, k, j)+ch3.at(i-1, k, jc)) @@ -503,7 +501,7 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } } else { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for k := 0; k < l1; k++ { for i := 2; i < ido; i += 2 { c13.set(i-1, k, j, ch3.at(i-1, k, j)+ch3.at(i-1, k, jc)) @@ -517,7 +515,7 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for k := 0; k < l1; k++ { c13.set(0, k, j, ch3.at(0, k, j)+ch3.at(0, k, jc)) c13.set(0, k, jc, ch3.at(0, k, jc)-ch3.at(0, k, j)) @@ -526,7 +524,7 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { ar1 := 1.0 ai1 := 0.0 for l := 1; l < ipph; l++ { - lc := ipp2 - (l + 1) + lc := ip - l ar1h := dcp*ar1 - dsp*ai1 ai1 = dcp*ai1 + dsp*ar1 ar1 = ar1h @@ -539,7 +537,7 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { ar2 := ar1 ai2 := ai1 for j := 2; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j ar2h := dc2*ar2 - ds2*ai2 ai2 = dc2*ai2 + ds2*ar2 ar2 = ar2h @@ -569,7 +567,7 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } } for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j j2 := 2 * j for k := 0; k < l1; k++ { cc3.set(ido-1, j2-1, k, ch3.at(0, k, j)) @@ -582,10 +580,10 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } if nbd < l1 { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j j2 := 2 * j for i := 2; i < ido; i += 2 { - ic := idp2 - (i + 1) + ic := ido - i for k := 0; k < l1; k++ { cc3.set(i-1, j2, k, ch3.at(i-1, k, j)+ch3.at(i-1, k, jc)) cc3.set(ic-1, j2-1, k, ch3.at(i-1, k, j)-ch3.at(i-1, k, jc)) @@ -597,7 +595,7 @@ func radfg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { return } for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j j2 := 2 * j for k := 0; k < l1; k++ { for i := 2; i < ido; i += 2 { @@ -764,8 +762,8 @@ func radb2(ido, l1 int, cc, ch, wa1 []float64) { } } for k := 0; k < l1; k++ { - ch3.set(ido-1, k, 0, cc3.at(ido-1, 0, k)+cc3.at(ido-1, 0, k)) - ch3.set(ido-1, k, 1, -(cc3.at(0, 1, k) + cc3.at(0, 1, k))) + ch3.set(ido-1, k, 0, 2*cc3.at(ido-1, 0, k)) + ch3.set(ido-1, k, 1, -2*cc3.at(0, 1, k)) } } @@ -963,8 +961,6 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { dcp := math.Cos(arg) dsp := math.Sin(arg) ipph := (ip + 1) / 2 - ipp2 := ip + 1 - idp2 := ido + 1 nbd := (ido - 1) / 2 if ido < l1 { @@ -982,7 +978,7 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j j2 := 2 * j for k := 0; k < l1; k++ { ch3.set(0, k, j, cc3.at(ido-1, j2-1, k)+cc3.at(ido-1, j2-1, k)) @@ -993,10 +989,10 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { if ido != 1 { if nbd < l1 { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j j2 := 2 * j for i := 2; i < ido; i += 2 { - ic := idp2 - (i + 1) + ic := ido - i for k := 0; k < l1; k++ { ch3.set(i-1, k, j, cc3.at(i-1, j2, k)+cc3.at(ic-1, j2-1, k)) ch3.set(i-1, k, jc, cc3.at(i-1, j2, k)-cc3.at(ic-1, j2-1, k)) @@ -1007,11 +1003,11 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } } else { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j j2 := 2 * j for k := 0; k < l1; k++ { for i := 2; i < ido; i += 2 { - ic := idp2 - (i + 1) + ic := ido - i ch3.set(i-1, k, j, cc3.at(i-1, j2, k)+cc3.at(ic-1, j2-1, k)) ch3.set(i-1, k, jc, cc3.at(i-1, j2, k)-cc3.at(ic-1, j2-1, k)) ch3.set(i, k, j, cc3.at(i, j2, k)-cc3.at(ic, j2-1, k)) @@ -1025,7 +1021,7 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { ar1 := 1.0 ai1 := 0.0 for l := 1; l < ipph; l++ { - lc := ipp2 - (l + 1) + lc := ip - l ar1h := dcp*ar1 - dsp*ai1 ai1 = dcp*ai1 + dsp*ar1 ar1 = ar1h @@ -1038,7 +1034,7 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { ar2 := ar1 ai2 := ai1 for j := 2; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j ar2h := dc2*ar2 - ds2*ai2 ai2 = dc2*ai2 + ds2*ar2 ar2 = ar2h @@ -1055,7 +1051,7 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } } for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for k := 0; k < l1; k++ { ch3.set(0, k, j, c13.at(0, k, j)-c13.at(0, k, jc)) ch3.set(0, k, jc, c13.at(0, k, j)+c13.at(0, k, jc)) @@ -1065,7 +1061,7 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { if ido != 1 { if nbd < l1 { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for i := 2; i < ido; i += 2 { for k := 0; k < l1; k++ { ch3.set(i-1, k, j, c13.at(i-1, k, j)-c13.at(i, k, jc)) @@ -1077,7 +1073,7 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { } } else { for j := 1; j < ipph; j++ { - jc := ipp2 - (j + 1) + jc := ip - j for k := 0; k < l1; k++ { for i := 2; i < ido; i += 2 { ch3.set(i-1, k, j, c13.at(i-1, k, j)-c13.at(i, k, jc)) @@ -1109,7 +1105,7 @@ func radbg(ido, ip, l1, idl1 int, cc, c1, c2, ch, ch2, wa []float64) { for k := 0; k < l1; k++ { idij := is for i := 2; i < ido; i += 2 { - idij = idij + 2 + idij += 2 c13.set(i-1, k, j, wa[idij-1]*ch3.at(i-1, k, j)-wa[idij]*ch3.at(i, k, j)) c13.set(i, k, j, wa[idij-1]*ch3.at(i, k, j)+wa[idij]*ch3.at(i-1, k, j)) }