| 1197 | } |
| 1198 | |
| 1199 | void cftmdl1(int n, float *a) |
| 1200 | { |
| 1201 | int i, i0, j, j0, j1, j2, j3, m, mh; |
| 1202 | float ew, w1r, w1i, wk1r, wk1i, wk3r, wk3i, |
| 1203 | wd1r, wd1i, wd3r, wd3i, ss1, ss3; |
| 1204 | float x0r, x0i, x1r, x1i, x2r, x2i, x3r, x3i; |
| 1205 | |
| 1206 | mh = n >> 3; |
| 1207 | m = 2 * mh; |
| 1208 | j1 = m; |
| 1209 | j2 = j1 + m; |
| 1210 | j3 = j2 + m; |
| 1211 | x0r = a[0] + a[j2]; |
| 1212 | x0i = a[1] + a[j2 + 1]; |
| 1213 | x1r = a[0] - a[j2]; |
| 1214 | x1i = a[1] - a[j2 + 1]; |
| 1215 | x2r = a[j1] + a[j3]; |
| 1216 | x2i = a[j1 + 1] + a[j3 + 1]; |
| 1217 | x3r = a[j1] - a[j3]; |
| 1218 | x3i = a[j1 + 1] - a[j3 + 1]; |
| 1219 | a[0] = x0r + x2r; |
| 1220 | a[1] = x0i + x2i; |
| 1221 | a[j1] = x0r - x2r; |
| 1222 | a[j1 + 1] = x0i - x2i; |
| 1223 | a[j2] = x1r - x3i; |
| 1224 | a[j2 + 1] = x1i + x3r; |
| 1225 | a[j3] = x1r + x3i; |
| 1226 | a[j3 + 1] = x1i - x3r; |
| 1227 | wd1r = 1; |
| 1228 | wd1i = 0; |
| 1229 | wd3r = 1; |
| 1230 | wd3i = 0; |
| 1231 | ew = M_PI_2 / m; |
| 1232 | w1r = (float)cos(2 * ew); |
| 1233 | w1i = (float)sin(2 * ew); |
| 1234 | wk1r = w1r; |
| 1235 | wk1i = w1i; |
| 1236 | ss1 = 2 * w1i; |
| 1237 | wk3i = 2 * ss1 * wk1r; |
| 1238 | wk3r = wk1r - wk3i * wk1i; |
| 1239 | wk3i = wk1i - wk3i * wk1r; |
| 1240 | ss3 = 2 * wk3i; |
| 1241 | i = 0; |
| 1242 | for (;;) { |
| 1243 | i0 = i + 4 * CDFT_LOOP_DIV; |
| 1244 | if (i0 > mh - 4) { |
| 1245 | i0 = mh - 4; |
| 1246 | } |
| 1247 | for (j = i + 2; j < i0; j += 4) { |
| 1248 | wd1r -= ss1 * wk1i; |
| 1249 | wd1i += ss1 * wk1r; |
| 1250 | wd3r -= ss3 * wk3i; |
| 1251 | wd3i += ss3 * wk3r; |
| 1252 | j1 = j + m; |
| 1253 | j2 = j1 + m; |
| 1254 | j3 = j2 + m; |
| 1255 | x0r = a[j] + a[j2]; |
| 1256 | x0i = a[j + 1] + a[j2 + 1]; |