Como tomo a FFT de dados com espaçamento desigual?

O algoritmo Fast Fourier Transform calcula uma decomposição de Fourier sob a suposição de que seus pontos de entrada estão igualmente espaçados no domínio do tempo, . E se não estiverem? Existe outro algoritmo que eu poderia usar, ou alguma forma de modificar a FFT, para explicar o que é efetivamente uma taxa de amostragem variável?$t_k = kT$

Se a solução depende de como as amostras são distribuídas, há duas situações específicas nas quais estou mais interessado:

• Taxa de amostragem constante com jitter: que é uma variável distribuída aleatoriamente. Suponha que seja seguro dizer .$t_k = kT + \delta t_k$$\delta t_k$$|\delta t_k| < T/2$
• Amostras descartadas de uma taxa de amostragem constante, de outra forma constante: que$t_k = n_k T$$n_k \in\mathbb{Z}\ge k$

Motivação: em primeiro lugar, essa foi uma das perguntas mais votadas na proposta deste site. Além disso, há algum tempo, participei de uma discussão sobre o uso da FFT (solicitada por uma pergunta no Stack Overflow ) na qual surgiram alguns dados de entrada com pontos de amostra desiguais. Aconteceu que os carimbos de data e hora estavam errados, mas isso me fez pensar em como alguém poderia resolver esse problema.

Existe uma grande variedade de técnicas para FFT não uniforme, e as mais eficientes são todas exatamente para o seu caso: amostras quase uniformes. A idéia básica é espalhar as fontes desigualmente amostradas em uma grade uniforme um pouco mais fina ("super-amostrada") através de convoluções locais contra gaussianos. Uma FFT padrão pode então ser executada na grade uniforme superamostrada e, em seguida, a convolução contra os gaussianos pode ser desfeita. Boas implementações são algo como vezes mais caro que uma FFT padrão em dimensões , onde é algo próximo de 4 ou 5.$C^d$$d$$C$

Eu recomendo ler Acelerando a transformação rápida de Fourier não uniforme de Greengard e Lee.

Também existem técnicas rápidas, isto é, ou mais rápidas, quando as fontes e / ou pontos de avaliação são escassos, e também há generalizações para operadores integrais mais gerais, por exemplo, Operadores Integrais de Fourier. Se você estiver interessado nessas técnicas, recomendo a transformação de Sparse Fourier via algoritmo de borboleta e um algoritmo rápido de borboleta para o cálculo de operadores integrais de Fourier . O preço pago nessas técnicas em relação às FFTs padrão é um coeficiente muito maior. Isenção de responsabilidade: meu orientador escreveu / escreveu esses dois trabalhos e passei bastante tempo paralelizando essas técnicas.$O(N^d \log N)$

Um ponto importante é que todas as técnicas acima são aproximações que podem ser arbitrariamente precisas às custas de tempos de execução mais longos, enquanto o algoritmo FFT padrão é exato.

No processamento do sinal, o aliasing é evitado enviando um sinal através de um filtro passa-baixo antes da amostragem. Jack Poulson já explicou uma técnica para FFT não uniforme usando gaussianos truncados como filtros de passa-baixo. Uma característica inconveniente dos gaussianos truncados é que, mesmo depois de ter decidido o espaçamento da grade para a FFT (= a taxa de amostragem no processamento do sinal), você ainda tem dois parâmetros livres: a largura do raio gaussiano e do truncamento.

Portanto, prefiro a função "hat" com uma largura de duas células da grade como filtro passa-baixo. Isso tem o efeito de que a ordem zero de Fourier é exata e que as ordens inferiores de Fourier convergirão quadraticamente. A transformação de Fourier da função "hat" é fácil de calcular (é o quadrado da função sinc), o que simplifica a destruição da convolução após a FFT. Observe que a função "hat" é a convolução da função característica da célula unitária (centralizada) consigo mesma. Qualquer taxa de convergência desejada pode ser alcançada, convoluindo a célula unitária mais de uma vez consigo mesma e usando a função resultante em vez da função "hat".

Se você estiver interessado em software, posso recomendar a biblioteca NFFT (em C com uma interface para o MATLAB), que pode ser encontrada aqui . Observe que também há uma biblioteca PFFT para computação paralela de FFT e até mesmo uma biblioteca PNFFT para FFTs paralelas não equispaces pelos mesmos desenvolvedores .

Além da resposta aceita. Aqui está um link para uma implementação de código aberto do método de Greengard e Lee: https://finufft.readthedocs.io/en/latest/ Possui wrappers para C, fortran, MATLAB, oitava e python. Eu acredito que o FINUFFT está escrito em C ++.

Ele é mantido e utilizado no Instituto NYU Courant, SFU, Instituto Flatiron (obviamente), Universidade do Texas Austin e Universidade Estadual da Flórida. Pelo menos estes são os que eu conheço.

Eu mesmo estou usando uma versão mais antiga, porque sou preguiçoso. Veja: https://cims.nyu.edu/cmcl/nufft/nufft.html

De interesse pode ser a Transformada de Fourier Discreta Compensada em Data:

Ferraz-Mello, S., 1981, Estimativa de períodos a partir de observações desigualmente espaçadas , The Astronomical Journal, 302: 757-763 .