Note that there are some explanatory texts on larger screens.

plurals
  1. PO
    primarykey
    data
    text
    <p>After the following subroutine</p> <pre><code>vDSP_ctoz((DSPComplex *)realData, 2, complexData, 1,samplesOver2); </code></pre> <p>is executed, <code>complexData</code> has <code>samplesOver2</code> elements. But soon after that, you call</p> <pre><code>vDSP_zvabs(complexData, 1, spectrumData, 1, samples); </code></pre> <p>which expects <code>complexData</code> to have <code>samples</code> elements, i.e. twice as many. This cannot be.</p> <p>Also, how is <code>realData</code> laid out? I ask because <code>vDSP_ctoz</code> expects its first argument to be laid out in the form</p> <pre><code>real0, imag0, real1, imag1, ... real(n-1), imag(n-1). </code></pre> <p>If your data is indeed real, then <code>imag0, imag1, ... imag(n-1)</code> should all be 0. If it is not, then <code>vDSP_ctoz</code> may not be expecting that. (Unless you are packing the real data in some clever way, which would be two [sic] clever by half!)</p> <p>Finally, <code>vDSP_create_fftsetup( 16, 2);</code> should probably be changed to</p> <pre><code>vDSP_create_fftsetup(16, 0); </code></pre> <p>===================================================================</p> <p>My sample code appended in postscript:</p> <pre><code> FFTSetup fftSetup = vDSP_create_fftsetup( 16, // vDSP_Length __vDSP_log2n, kFFTRadix2 // FFTRadix __vDSP_radix // CAUTION: kFFTRadix2 is an enum that is equal to 0 // kFFTRadix5 is an enum that is equal to 2 // DO NOT USE 2 IF YOU MEAN kFFTRadix2 ); NSAssert(fftSetup != NULL, @"vDSP_create_fftsetup() failed to allocate storage"); int numSamples = 65536; // numSamples must be an integer power of 2; in this case 65536 = 2 ^ 16 float realData[numSamples]; // Prepare the real data with (ahem) fake data, in this case // the sum of 3 sinusoidal waves representing a C major chord. // The fake data is rigged to have a sampling frequency of 44100 Hz (as for a CD). // As always, the Nyquist frequency is just half the sampling frequency, i.e., 22050 Hz. for (int i = 0; i &lt; numSamples; i++) { realData[i] = sin(2 * M_PI * 261.76300048828125 * i / 44100.0) // C4 = 261.626 Hz + sin(2 * M_PI * 329.72717285156250 * i / 44100.0) // E4 = 329.628 Hz + sin(2 * M_PI * 392.30804443359375 * i / 44100.0); // G4 = 391.995 Hz } float splitReal[numSamples / 2]; float splitImag[numSamples / 2]; DSPSplitComplex splitComplex; splitComplex.realp = splitReal; splitComplex.imagp = splitImag; vDSP_ctoz( (const DSPComplex *)realData, // const DSPComplex __vDSP_C[], 2, // vDSP_Stride __vDSP_strideC, MUST BE A MULTIPLE OF 2 &amp;splitComplex, // DSPSplitComplex *__vDSP_Z, 1, // vDSP_Stride __vDSP_strideZ, (numSamples / 2) // vDSP_Length __vDSP_size ); vDSP_fft_zrip( fftSetup, // FFTSetup __vDSP_setup, &amp;splitComplex, // DSPSplitComplex *__vDSP_ioData, 1, // vDSP_Stride __vDSP_stride, (vDSP_Length)lround(log2(numSamples)), // vDSP_Length __vDSP_log2n, // IMPORTANT: THE PRECEDING ARGUMENT MUST BE LOG_BASE_2 OF THE NUMBER OF floats IN splitComplex // FOR OUR EXAMPLE, THIS WOULD BE (numSamples / 2) + (numSamples / 2) = numSamples kFFTDirection_Forward // FFTDirection __vDSP_direction ); printf("DC component = %f\n", splitComplex.realp[0]); printf("Nyquist component = %f\n\n", splitComplex.imagp[0]); // Next, we compute the Power Spectral Density (PSD) from the FFT. // (The PSD is just the magnitude-squared of the FFT.) // (We don't bother with scaling as we are only interested in relative values of the PSD.) float powerSpectralDensity[(numSamples / 2) + 1]; // the "+ 1" is to make room for the Nyquist component // We move the Nyquist component out of splitComplex.imagp[0] and place it // at the end of the array powerSpectralDensity, squaring it as we go: powerSpectralDensity[numSamples / 2] = splitComplex.imagp[0] * splitComplex.imagp[0]; // We can now zero out splitComplex.imagp[0] since the imaginary part of the DC component is, in fact, zero: splitComplex.imagp[0] = 0.0; // Finally, we compute the squares of the magnitudes of the elements of the FFT: vDSP_zvmags( &amp;splitComplex, // DSPSplitComplex *__vDSP_A, 1, // vDSP_Stride __vDSP_I, powerSpectralDensity, // float *__vDSP_C, 1, // vDSP_Stride __vDSP_K, (numSamples / 2) // vDSP_Length __vDSP_N ); // We print out a table of the PSD as a function of frequency // Replace the "&lt; 600" in the for-loop below with "&lt;= (numSamples / 2)" if you want // the entire spectrum up to and including the Nyquist frequency: printf("Frequency_in_Hz Power_Spectral_Density\n"); for (int i = 0; i &lt; 600; i++) { printf("%f, %f\n", (i / (float)(numSamples / 2)) * 22050.0, powerSpectralDensity[i]); // Recall that the array index i = 0 corresponds to zero frequency // and that i = (numSamples / 2) corresponds to the Nyquist frequency of 22050 Hz. // Frequency values intermediate between these two limits are scaled proportionally (linearly). } // The output PSD should be zero everywhere except at the three frequencies // corresponding to the C major triad. It should be something like this: /*************************************************************************** DC component = -0.000000 Nyquist component = -0.000000 Frequency_in_Hz Power_Spectral_Density 0.000000, 0.000000 0.672913, 0.000000 1.345825, 0.000000 2.018738, 0.000000 2.691650, 0.000000 . . . 260.417175, 0.000000 261.090088, 0.000000 261.763000, 4294967296.000000 262.435913, 0.000000 263.108826, 0.000000 . . . 328.381348, 0.000000 329.054260, 0.000000 329.727173, 4294967296.000000 330.400085, 0.000000 331.072998, 0.000000 . . . 390.962219, 0.000000 391.635132, 0.000000 392.308044, 4294966784.000000 392.980957, 0.000000 393.653870, 0.000000 . . . ***************************************************************************/ vDSP_destroy_fftsetup(fftSetup); </code></pre>
    singulars
    1. This table or related slice is empty.
    plurals
    1. This table or related slice is empty.
    1. This table or related slice is empty.
    1. This table or related slice is empty.
    1. VO
      singulars
      1. This table or related slice is empty.
    2. VO
      singulars
      1. This table or related slice is empty.
    3. VO
      singulars
      1. This table or related slice is empty.
 

Querying!

 
Guidance

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload