Ordem das frequências do MATLAB FFT

9

Este wikibook afirma que a saída do MATLAB FFTcorresponde aos números de onda ordenados como:

k={0,1,...,n2,n2+1,n2+2,...,1}

No entanto, nos códigos de exemplo da mesma página, os números de onda são codificados como

k = [0:n/2-1 0 -n/2+1:-1];

que é o mesmo que o primeiro, mas com o número onda (o "número máximo de ondas") substituído por . Parece estranho que eles incluam duas vezes.0 0n/200

Parece que a ordem correta é necessária para obter derivadas através de transformadas de Fourier, conforme descrito no wikibook. Qual delas está correta e o MATLAB documenta isso em qualquer lugar?

Dúvida
fonte
3
veja a função fftshift. Ele pegará a saída de fft e a reordenará de [-n / 2 + 1: n / 2-1], o que deve ajudar com sua confusão.
Godric Seer
Isso não é algo que você pode testar rapidamente? Descobrir em que ordem a saída da FFT é fácil de determinar com algumas experiências.
Federico Poloni

Respostas:

4

Desejo expandir meu comentário e refazer o exemplo que você menciona de uma maneira que seja mais compreensível que o original e explicar por que fftretorna os coeficientes da maneira que é.

Para referência, a parte fft do exemplo é:

Nx = size(x,2);
k = 2*pi/(b-a)*[0:Nx/2-1 0 -Nx/2+1:-1];
dFdx = ifft(1i*k.*fft(f));
d2Fdx2 = ifft(-k.^2.*fft(f));

Adicionei outra seção de código diretamente abaixo:

Nx = size(x,2);
k = 2*pi/(b-a)*(-Nx/2:Nx/2-1);
dFdxp = ifft(ifftshift(1i*k.*fftshift(fft(f))));
d2Fdx2p = ifft(ifftshift(-k.^2.*fftshift(fft(f))));

e embrulhou os dois pedaços de código em um período tic; tocde tempo necessário. Em um formato mais legível, o segundo método usa:

ckf = fftshift(fft(f));
ckdf = 1i*k.*ckf;
df = ifft(ifftshift(ckdf));

A primeira diferença é que o segundo exemplo é muito mais intuitivo k. Essa é a principal vantagem do segundo exemplo, já que k está agora na forma que pensamos sobre eles. Nas segunda e terceira linhas, eu tinha que adicionar fftshifta chamada para ffte depois uma chamada ifftshiftdiretamente para dentro da chamada ifft. Essas chamadas de função adicionais reordenam os coeficientes do que é necessário para o computador trabalhar com eles da maneira que os humanos geralmente pensam sobre eles.

A questão do segundo exemplo é que, embora kseja mais intuitivo para nós, isso deixa as matrizes internas para resolver e inverter fftem formas que não são tão vantajosas. Portanto, ou temos que mudar a ordem com chamadas para fftswitche ifftswitchou deve ser codificado para as fftfunções. Isso é menos propenso a erros por parte dos usuários (supondo que eles não estejam familiarizados com o funcionamento do fft, como muitas pessoas o fazem), mas você paga um preço em tempo de execução.

Como afirmei antes, adicionei chamadas de tempo nos dois blocos para comparação e corri para vários N. Os resultados de tempo foram:

N =     1000,  Ex1 = 0.000222 s,   Ex2 = 0.007072 s
N =    10000,  Ex1 = 0.001576 s,   Ex2 = 0.003506 s
N =   100000,  Ex1 = 0.023857 s,   Ex2 = 0.034051 s
N =  1000000,  Ex1 = 0.213816 s,   Ex2 = 0.406250 s
N = 10000000,  Ex1 = 4.555143 s,   Ex2 = 7.102348 s

Como você pode ver, o ato de alternar os valores diminui consideravelmente o processo, especialmente em N baixo (onde é 30x mais lento). Este é apenas um exemplo, e seu computador pode mostrar tendências ligeiramente diferentes, dependendo de coisas como velocidade da memória, núcleos / velocidade do processador, etc., mas é ilustrativo do ponto. A razão pela qual a fftsaída é confusa é porque está economizando uma fração não trivial do seu tempo de computação.

Godric Seer
fonte
3
Isso ainda não responde à pergunta original - por que funciona para incluir zero duas vezes? E por que a documentação sugere que o zero deve ser incluído apenas uma vez?
David Ketcheson
2

Sua pergunta sobre a substituição do número de onda é bastante complicada. Em geral, essa modificação do número de onda não tem como objetivo salvar os fracassos, como alguns sugeriram aqui, mas projetada para respeitar as peculiaridades analíticas de, digamos, certos operadores diferenciais. Estou muito surpreso por não ter encontrado uma discussão relacionada nos Métodos espectrais de Trefethen. Para prosseguir, assumirei que você está preocupado com a diferenciação espectral baseada em FFT e que está realizando transformações em um domínio de cardinalidade uniforme.

A regra geral, para derivadas ímpares, é definir

k = [0:n/2-1 0 -n/2+1:-1];

e mesmo para derivativos, para definir

k = [0:n/2 -n/2+1:-1];

Se você estiver lidando com algum preenchimento, o tratamento da frequência de Nyquist é igualmente tedioso. Há um excelente artigo que comprova e observa esses tópicos aqui !

Stuart Hilton
fonte
1

n/2nn/20

Para uma breve explicação e demonstração do alias, consulte o meu caderno IPython aqui .

David Ketcheson
fonte
0

NN/2xRN/2N/2N/2

mu7z
fonte