Function I wrote that optimizes the DFT or iDFT for BOTH prime and composite signals...
The forward transform is triggered by -1 and takes the time-signal as a row vector:
[Y] = fftmodule(y,-1)
The inverse transform auto-normalizes by N, is triggered by 1 and takes the frequency-signal as a row vector:
[y] = fftmodule(Y,1)
The algorithm decimates to N's prime factorization following the branches and nodes of a factor tree. In formal literature this may be referred to as Mixed Radix FFT, but its really just recursive decimation of additive groups and this method is easily derivable via circular convolutions :)
At the prime tree level, algorithm either performs a naive DFT or if needed performs a single Rader's Algorithm Decomposition to (M-1), zero-pads to power-of-2, then proceeds to Rader's Convolution routine (not easy to derive). Finally it upsamples through the origianlly strucutured nodes and branches incorporating twiddle factors for the solution.