Duck
Duck

Reputation: 35933

Doing a Discrete Fourier Transformed followed by an Inverse Discrete Fourier Transform does not produce the original values on Accelerate vDSP

I have a series of discrete values. So, I have no imaginary values for the input part. I am doing a Discrete Fourier Transform to these values and performing an Inverse Discrete Fourier Transform, to get the same values back. The idea was to test if vDSP from Accelerate Framework was working as expected, cause the documentation is a piece of crap and you have to figure it out by yourself.

Consider the following code

  lazy var DFTSetupForward: vDSP_DFT_Setup = {
    guard let setup = vDSP_DFT_zop_CreateSetupD(
      nil,
      vDSP_Length(self.numSamples),
      vDSP_DFT_Direction.FORWARD) else {
        fatalError("can't create vDSP_DFT_Setup")
    }
    return setup
  }()

  lazy var DFTSetupInverse: vDSP_DFT_Setup = {
    guard let setup = vDSP_DFT_zop_CreateSetupD(
      nil,
      vDSP_Length(self.numSamples),
      vDSP_DFT_Direction.INVERSE) else {
        fatalError("can't create vDSP_DFT_Setup")
    }
    return setup
  }()
 func discreteFourierTransform (_ valores:[Double] = []) -> ([Double], [Double]) {
    let numeroDados = valores.count

    let inputImag  = Array<Double>(repeating:0.0, count:numeroDados)
    var outputReal = Array<Double>(repeating:0.0, count:numeroDados)
    var outputImag = Array<Double>(repeating:0.0, count:numeroDados)

    vDSP_DFT_ExecuteD(DFTSetupForward, valores, inputImag, &outputReal, &outputImag)
    return (outputReal, outputImag)
  }


  // faz a DISCRETE FOURIER TRANSFORM DOUBLE
  // a saída corresponde aos vetores reais e imaginários
  func discreteFourierTransformInverse (_ valoresReais:[Double] = [],
                                        _ valoresImaginarios:[Double] = []) -> ([Double], [Double]) {
    let numeroDados = valoresReais.count

    var outputReal = Array<Double>(repeating:0.0, count:numeroDados)
    var outputImag = Array<Double>(repeating:0.0, count:numeroDados)

    vDSP_DFT_ExecuteD(DFTSetupInverse, valoresReais, valoresImaginarios, &outputReal, &outputImag)
    return (outputReal, outputImag)
  }

and later these calls

let (outputDFTreal, outputDFTImaginario ) = self.discreteFourierTransform(normalizedData)

let (sinalRealX, sinalImaginarioX ) = self.discreteFourierTransformInverse(outputDFTreal, outputDFTImaginario)

or in other words, I am taking the results of the Fourier Transform from the first let and injecting them into the Inverse Fourier Transform. I expect sinalRealX to be equal to the original data, that in this case is normalizedData and also expect signalImaginarioX to be all zeros.

But the values are completely different!

What is going on?

Upvotes: 2

Views: 578

Answers (2)

hotpaw2
hotpaw2

Reputation: 70673

Floating point numbers are an approximation. And 10^-19 is approximately zero, relative to cos(0), in double precision floating-point arithmetic.

So you did get all zeros in the result.

Upvotes: 2

Rob Napier
Rob Napier

Reputation: 299265

You shouldn't expect the IDFT to return the normalized values. You should expect it to return the orginal (pre-normalize) values, which I believe it does (within computational error):

let fft = FFT()
let data = Array(stride(from: 0.0, to: Double(fft.numSamples), by: 1))
let normalizedData = data.map{ $0/Double(fft.numSamples) }
let (outputDFTreal, outputDFTImaginario ) = fft.discreteFourierTransform(normalizedData)
let (sinalRealX, sinalImaginarioX ) = fft.discreteFourierTransformInverse(outputDFTreal, outputDFTImaginario)
print(data)
print(sinalRealX)
print(sinalImaginarioX)

let zeros = Array(repeating: 0.0, count: data.count)
func avgerr(_ a: [Double], _ b: [Double]) -> Double {
    return zip(a, b).map { abs($0 - $1) }.reduce(0, +) / Double(a.count)
}

avgerr(data, sinalRealX)  // 2.442490654175344e-15
avgerr(zeros, sinalImaginarioX) // 2.602442788260593e-15

Upvotes: 3

Related Questions