SHOW:
|
|
- or go back to the newest paste.
| 1 | # modulo fft | |
| 2 | from math import * | |
| 3 | ||
| 4 | def fft(m, f=[]): #genera fft | |
| 5 | result = [] | |
| 6 | - | n=2<<(m-1) |
| 6 | + | n=2<<(m-1) |
| 7 | nesp=0.5 | |
| 8 | for k in range(1,m+1): | |
| 9 | d=n/(2<<(k-1) | |
| 10 | nesp=nesp*2 | |
| 11 | s=0 | |
| 12 | for i in range(1,nesp+1): | |
| 13 | r=s | |
| 14 | z1=r>>(m-k) | |
| 15 | z=0 | |
| 16 | for q in range(k): | |
| 17 | z=(((z1>>q) & 01) | (z<<1) | |
| 18 | ||
| 19 | w=2*pi*z/(2<<(k-1)) | |
| 20 | e=cos(w)-sin(w)j | |
| 21 | for l in range (d): | |
| 22 | r=l+s | |
| 23 | r1=r+d | |
| 24 | t=e*f[r1] | |
| 25 | f[r1]=f[r]-t | |
| 26 | f[r]=f[r]+t | |
| 27 | s=s+2*d | |
| 28 | #bitReversal(n, m, f=[]) | |
| 29 | for x in range(1, n-1): | |
| 30 | y=0 | |
| 31 | for i in range(m): | |
| 32 | y=(((x>>i) & 01) | (y<<1)) | |
| 33 | if x<y : | |
| 34 | t=f[x] | |
| 35 | f[x]=f[y] | |
| 36 | f[y]=t | |
| 37 | #Ampiezze fasi | |
| 38 | nmezzi=n/2 | |
| 39 | for i in range(n): | |
| 40 | f[i]=f[i]/nmezzi | |
| 41 | A[0]=f[0].real | |
| 42 | phi[0]=0 | |
| 43 | for k in range(1, nmezzi+1): | |
| 44 | A[k]=sqrt(f[k].real*f[k].real+f[k].imag*f[k].imag]) | |
| 45 | result.append(A[k]) | |
| 46 | if f[k].real == 0: | |
| 47 | if f[k].imag>0: | |
| 48 | phi[k]=pi | |
| 49 | if f[k].imag==0: | |
| 50 | phi[k]=0 | |
| 51 | if f[k].imag<0: | |
| 52 | phi[k]=2*pi-pi/2 | |
| 53 | if f[k].real>0: | |
| 54 | if f[k]>=0: | |
| 55 | phi[k]=atan(f[k].imag, f[k].real) | |
| 56 | if f[k]<0: | |
| 57 | phi[k]=2*pi + atan(f[k].imag, f[k].real) | |
| 58 | if f[k].real<0: | |
| 59 | phi[k]=atan(f[k].imag, f[k].real) | |
| 60 | return result |
