Utiliser la pomme FFT et accélérer le framework

Question

Quelqu'un a-t-il déjà utilisé le apple FFT pour une application iPhone ou sait-il où je pourrais trouver un exemple d'application pour savoir comment l'utiliser? Je sais que Apple a affiché un exemple de code, mais je ne sais pas trop comment le mettre en œuvre dans un projet réel.

Answer 1

Je viens de recevoir la FFT code de travail pour un iPhone de projet:

• créer un nouveau projet
• supprimez tous les fichiers sauf pour le principal.m et xxx_info.plist
• allez à paramètres du projet et de recherche pour la pch et de l'empêcher d'essayer de charger une .pch (comme nous venons de le supprimé)
• copier coller le code exemple de plus de ce que vous avez en main.m
• supprimer la ligne #include Carbone. Le carbone est OSX.
• supprimer tous les cadres, et ajouter accélérer le cadre

Vous pourriez aussi avoir besoin de supprimer une entrée de info.plist que dit le projet de charger un xib, mais je suis sûr à 90% que vous n'avez pas besoin de s'embêter avec ça.

REMARQUE: le Programme des sorties de la console, les résultats sont comme 0.000 ce n'est pas une erreur, c'est juste très très vite

Ce code est vraiment bêtement obscur; il est généreusement commenté, mais les commentaires ne sont pas réellement rendre la vie plus facile.

Fondamentalement, le cœur du problème est:

vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_FORWARD);
vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_INVERSE);

FFT sur n réel des flotteurs, et puis en sens inverse pour revenir là où nous avons commencé. ip signifie en place, ce qui signifie &A est écrasée C'est la raison de tout cet emballage spécial sottises, de sorte que nous pouvons écraser la valeur de retour dans le même espace que de l'envoyer en valeur.

Pour donner une idée (comme: pourquoi serions-nous à l'aide de cette fonction, en premier lieu?), Disons que nous voulons effectuer la détection de hauteur sur entrée microphone, et nous avons le configurer de sorte que certains de rappel se déclenche à chaque fois que le microphone est en 1024 flotteurs. En supposant que le microphone de fréquence d'échantillonnage est de 44,1 kHz, de sorte que ~44 images / sec.

Donc, de notre temps-fenêtre quelle que soit la durée du temps de 1024 échantillons, c'est à dire 1/44 s.

Donc on aurait Un pack avec 1024 flotte à partir du microphone, réglez log2n=10 (2^10=1024), précalculer quelques bobines (setupReal) et:

vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_FORWARD);

Maintenant contiendra n/2 nombres complexes. Ceux-ci représentent les n/2 bacs de fréquences:

• bin[1].idealFreq = 44Hz -- c'est à dire La fréquence la plus basse que nous pouvons détecter de manière fiable est UNE onde à l'intérieur de cette fenêtre, c'est à dire un 44Hz vague.

• bin[2].idealFreq = 2 * 44Hz

• etc.

• bin[512].idealFreq = 512 * 44Hz -- la fréquence La plus élevée, nous pouvons détecter (connu comme la fréquence de Nyquist) est l'endroit où chaque paire de points représente une onde, c'est à dire 512 complète des ondes à l'intérieur de la fenêtre, c'est à dire 512 * 44Hz, ou: n/2 * bin[1].idealFreq

• En fait il y a un Bin, Bin[0], qui est souvent désigné comme "DC Offset'. Il se trouve que Bin[0] et Ben[n/2] aura toujours complexe de la composante 0, A[0].realp est utilisé pour stocker Bin[0] et[0].imagp est utilisé pour stocker Bin[n/2]

Et l'importance de chaque nombre complexe est la quantité d'énergie de vibration autour de cette fréquence.

Donc, comme vous pouvez le voir, il ne serait pas un très grand terrain de détecteur car il n'a pas presque assez fine granularité. Il y a une ruse de l'Extraction précise des fréquences de FFT des Bacs à l'aide de changement de phase entre les images pour obtenir la fréquence exacte pour un bin.

Ok, Maintenant sur le code:

Notez que le 'ip' dans vDSP_fft_zrip, = 'au lieu de' ie de sortie remplace Un ("r" signifie qu'il prend un réel entrées)

Regardez la documentation sur vDSP_fft_zrip,

Réel les données sont stockées dans split complexe forme, avec de curieuses réels stockés sur le imaginaire côté de la division complexe forme et même de reals en stockées sur le côté réel.

c'est probablement la chose la plus difficile à comprendre. Nous utilisons le même conteneur (&A) tout le chemin à travers le processus. donc, au début, nous voulons remplir avec de n nombres réels. après la FFT, il va être tenue n/2 nombres complexes. nous avons ensuite jeter que dans la transformation inverse, et nous espérons sortir notre origine n nombres réels.

maintenant, la structure de son programme d'installation pour des valeurs complexes. Donc vDSP besoins de standardiser la façon d'emballer des nombres réels.

donc, nous avons d'abord générer des n nombres réels: 1, 2, ..., n

for (i = 0; i < n; i++)
    originalReal[i] = (float) (i + 1);

Ensuite, nous les emballer dans Un que n/2 complexe #s:

// 1. masquerades n real #s as n/2 complex #s = {1+2i, 3+4i, ...}
// 2. splits to 
//   A.realP = {1,3,...} (n/2 elts)
//   A.compP = {2,4,...} (n/2 elts)
//
vDSP_ctoz(
          (COMPLEX *) originalReal, 
          2,                            // stride 2, as each complex # is 2 floats
          &A, 
          1,                            // stride 1 in A.realP & .compP
          nOver2);                      // n/2 elts

Vous avez vraiment besoin de regarder la façon dont Un est alloué pour l'obtenir, peut-être rechercher COMPLEX_SPLIT dans la documentation.

A.realp = (float *) malloc(nOver2 * sizeof(float));
A.imagp = (float *) malloc(nOver2 * sizeof(float));

Ensuite, nous faire un pré-calcul.

---

Rapide DSP de la classe de mathématiques bods: Fourier théorie prend beaucoup de temps pour obtenir autour de votre tête (j'ai été regarder sur et en dehors depuis plusieurs années maintenant)

Un cisoid est:

z = exp(i.theta) = cos(theta) + i.sin(theta)

c'est à dire un point sur le cercle unité dans le plan complexe.

Quand vous multipliez les nombres complexes, les angles ajouter. Donc z^k va garder se déplaçant sur le cercle unité; z^k peuvent être trouvés à un angle de k.theta

• Choisissez z1 = 0+1i, c'est à dire un quart de tour de l'axe réel, et de l'avis que z1^2 z1^3 z1^4 donnent chacun un quart de tour de sorte que z1^4 = 1

• Choisissez z2 = -1, c'est à dire un demi-tour. aussi z2^4 = 1 et z2 a terminé 2 cycles à ce point (z2^2 est également = 1). Alors que vous pourriez penser z1 comme la fréquence fondamentale et z2 premier harmonique

• De même, z3 = les trois-quarts de la révolution " point c'est à dire -je termine exactement 3 cycles, mais en réalité, le futur 3/4 chaque fois c'est la même chose que d'aller en arrière 1/4 à chaque fois

c'est à dire z3 est juste z1, mais dans la direction opposée, qui s'appelle l'aliasing

z2 est le plus significatif de la fréquence, nous avons choisi 4 échantillons de tenir une vague complète.

• z0 = 1+0i, z0^(n'importe quoi)=1, c'est l'offset DC

Vous pouvez exprimer toutes les 4 points du signal comme une combinaison linéaire de z0, z1 et z2 c'est à dire que Vous êtes les projette sur ces vecteurs de base

mais je vous entends demander "que signifie pour le projet d'un signal sur un cisoid?"

Vous pouvez pensez-y de cette façon: L'aiguille tourne autour de la cisoid, donc, à l'exemple de k, l'aiguille est orientée dans la direction k.thêta, et la longueur du signal[k]. Un signal qui correspond à la fréquence de la cisoid exactement renflement la forme obtenue dans une certaine direction. Donc, si vous additionnez toutes les contributions, vous aurez un fort vecteur résultant. Si la fréquence est pratiquement à la hauteur, que le renflement sera plus faible et se déplacera lentement autour du cercle. Pour un signal qui ne correspond pas à la fréquence, les contributions annuler l'un l'autre.

http://complextoreal.com/tutorials/tutorial-4-fourier-analysis-made-easy-part-1/ vous aidera à obtenir une compréhension intuitive.

Mais l'essentiel est; si nous avons choisi de projet 1024 échantillons sur {z0,...,z512} nous aurions précalculer z0 thru z512, et c'est ce que ce calcul prévisionnel étape est.

---

Notez que si vous faites cela dans le code réel, vous voulez probablement faire cela qu'une seule fois lorsque l'application des charges et de l'appel de la complémentarité de la libération de la fonction une fois quand il se ferme. Ne le faites PAS beaucoup de temps -- c'est cher.

// let's say log2n = 8, so n=2^8=256 samples, or 'harmonics' or 'terms'
// if we pre-calculate the 256th roots of unity (of which there are 256) 
// that will save us time later.
//
// Note that this call creates an array which will need to be released 
// later to avoid leaking
setupReal = vDSP_create_fftsetup(log2n, FFT_RADIX2);

Il est intéressant de noter que si l'on log2n par exemple à 8, vous pouvez jeter ces valeurs précalculées dans toute la fft fonction qui utilise la résolution <= 2^8. Donc (sauf si vous voulez ultime de la mémoire optimisation) il suffit de créer un jeu pour les plus de haute résolution, vous allez avoir besoin, et de l'utiliser pour tout.

Maintenant, le réel se transforme, en faisant l'utilisation des choses que nous venons de précalculées:

vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_FORWARD);

À ce stade, Un contient n/2 nombres complexes, seul le premier est en fait deux nombres réels (l'offset DC Nyquist #) se faisant passer pour un nombre complexe. La présentation de la documentation explique cet emballage. Il est tout à fait soigné -- fondamentalement, il permet l' (complexe) des résultats de la transformation pour être emballé dans le même espace de mémoire que le (vrai, mais bizarrement) emballage d'entrées.

vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_INVERSE);

et de retour... nous aurons encore besoin de décompresser de notre tableau original de A. ensuite, nous comparons juste pour vérifier que nous avons obtenu exactement ce que nous avons commencé avec la libération de notre prédéterminée de bobines et fait!

Mais attendez! avant de déballer, il y a une dernière chose qui doit être fait:

// Need to see the documentation for this one...
// in order to optimise, different routines return values 
// that need to be scaled by different amounts in order to 
// be correct as per the math
// In this case...
scale = (float) 1.0 / (2 * n);

vDSP_vsmul(A.realp, 1, &scale, A.realp, 1, nOver2);
vDSP_vsmul(A.imagp, 1, &scale, A.imagp, 1, nOver2);

Answer 2

Voici un exemple concret: Un extrait de c ++ qui utilise les routines fft vDSP d’Accelerate pour effectuer une autocorrélation sur l’entrée de l’unité audio Remote IO. Utiliser ce framework est assez compliqué, mais la documentation n’est pas si mauvaise.

 OSStatus DSPCore::initialize (double _sampleRate, uint16_t _bufferSize) {
    sampleRate = _sampleRate;
    bufferSize = _bufferSize;
    peakIndex = 0;
    frequency = 0.f;
    uint32_t maxFrames = getMaxFramesPerSlice();
    displayData = (float*)malloc(maxFrames*sizeof(float));
    bzero(displayData, maxFrames*sizeof(float));
    log2n = log2f(maxFrames);
    n = 1 << log2n;
    assert(n == maxFrames);
    nOver2 = maxFrames/2;
    A.realp = (float*)malloc(nOver2 * sizeof(float));
    A.imagp = (float*)malloc(nOver2 * sizeof(float));
    FFTSetup fftSetup = vDSP_create_fftsetup(log2n, FFT_RADIX2);

    return noErr;
}

void DSPCore::Render(uint32_t numFrames, AudioBufferList *ioData) {

    bufferSize = numFrames;
    float ln = log2f(numFrames);

    //vDSP autocorrelation

    //convert real input to even-odd
    vDSP_ctoz((COMPLEX*)ioData->mBuffers[0].mData, 2, &A, 1, numFrames/2);
    memset(ioData->mBuffers[0].mData, 0, ioData->mBuffers[0].mDataByteSize);
    //fft
    vDSP_fft_zrip(fftSetup, &A, 1, ln, FFT_FORWARD);

    // Absolute square (equivalent to mag^2)
    vDSP_zvmags(&A, 1, A.realp, 1, numFrames/2);
    bzero(A.imagp, (numFrames/2) * sizeof(float));    

    // Inverse FFT
    vDSP_fft_zrip(fftSetup, &A, 1, ln, FFT_INVERSE);

    //convert complex split to real
    vDSP_ztoc(&A, 1, (COMPLEX*)displayData, 2, numFrames/2);

    // Normalize
    float scale = 1.f/displayData[0];
    vDSP_vsmul(displayData, 1, &scale, displayData, 1, numFrames);

    // Naive peak-pick: find the first local maximum
    peakIndex = 0;
    for (size_t ii=1; ii < numFrames-1; ++ii) {
        if ((displayData[ii] > displayData[ii-1]) && (displayData[ii] > displayData[ii+1])) {
            peakIndex = ii;
            break;
        }
    }

    // Calculate frequency
    frequency = sampleRate / peakIndex + quadInterpolate(&displayData[peakIndex-1]);

    bufferSize = numFrames;

    for (int ii=0; ii<ioData->mNumberBuffers; ++ii) {
        bzero(ioData->mBuffers[ii].mData, ioData->mBuffers[ii].mDataByteSize);
    }
}
 

Answer 3

Alors je vais dire Apple FFT Cadre est rapide... Vous avez besoin de savoir comment une FFT afin d'optimiser la détection de hauteur (c'est à dire le calcul de la différence de phase sur chaque FFT afin de trouver la hauteur exacte, pas la hauteur de la plus dominer bin).

Je ne sais pas si c'est de l'aide, mais j'ai téléchargé mon Terrain de Détecteur d'objet à partir de mon tuner app (musicianskit.com/developer.php). Il est un exemple xCode 4 projet de télécharger aussi (donc vous pouvez voir comment la mise en œuvre de travaux).

Je suis en train de travailler sur la mise en ligne d'un exemple de la FFT de la mise en œuvre -- alors restez à l'écoute et je vais mettre à jour ce une fois que cela arrive.

Amusez-vous bien!

Answer 4

Voici un autre exemple concret : https://github.com/krafter/DetectingAudioFrequency.