誰かがApple FFTiPhone アプリ用に使用したことがありますか、それとも使用方法に関するサンプル アプリケーションがどこにあるか知っていますか? Apple がいくつかのサンプル コードを投稿していることは知っていますが、それを実際のプロジェクトに実装する方法がよくわかりません。
55653 次
4 に答える
139
iPhone プロジェクトで動作する FFT コードを取得しました。
- 新しいプロジェクトを作成する
- main.m と xxx_info.plist を除くすべてのファイルを削除します
- プロジェクト設定に移動して pch を検索し、.pch のロードを停止します (削除したばかりであることを確認してください)。
- main.m にあるものにコード例をコピーして貼り付けます
- #include の Carbon の行を削除します。Carbon は OSX 用です。
- すべてのフレームワークを削除し、accelerated フレームワークを追加します
プロジェクトにxibをロードするように指示するinfo.plistからエントリを削除する必要があるかもしれませんが、それを気にする必要はないと90%確信しています。
注: プログラムはコンソールに出力します。結果は 0.000 として出力されますが、これはエラーではありません。非常に高速です。
このコードは本当にばかばかしいほどあいまいです。寛大にコメントされていますが、コメントは実際には生活を楽にするものではありません。
基本的にその中心は次のとおりです。
vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_FORWARD);
vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_INVERSE);
n 個の実数フロートで FFT を実行し、逆方向に実行して最初の場所に戻ります。ip は in-place の略で、これは &A が上書きされることを意味します。これが、このすべての特別なパッキング マラーキーの理由です。つまり、戻り値を送信値と同じスペースに押しつぶすことができます。
いくつかの視点を与えるために (たとえば、そもそもなぜこの関数を使用するのでしょうか?)、マイク入力でピッチ検出を実行したいとしましょう。毎回コールバックがトリガーされるように設定しました。マイクは 1024 浮動小数点で取得します。マイクのサンプリング レートが 44.1kHz だとすると、1 秒あたり約 44 フレームになります。
したがって、タイム ウィンドウは 1024 サンプルの持続時間、つまり 1/44 秒です。
したがって、マイクからの 1024 個のフロートを A にパックし、log2n=10 (2^10=1024) を設定し、いくつかのボビンを事前に計算し (setupReal)、次のようにします。
vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_FORWARD);
A には n/2 の複素数が含まれます。これらは n/2 周波数ビンを表します。
bin[1].idealFreq = 44Hz -- つまり、確実に検出できる最低周波数は、そのウィンドウ内の 1 つの完全な波、つまり 44Hz の波です。
bin[2].idealFreq = 2 * 44Hz
等
bin[512].idealFreq = 512 * 44Hz -- 検出できる最高周波数 (ナイキスト周波数と呼ばれます) は、ポイントのすべてのペアが波を表す場所です。つまり、ウィンドウ内の 512 個の完全な波、つまり 512 * 44Hz、または: n/2 * bin[1].idealFreq
実際には、「DC オフセット」と呼ばれることが多い Bin[0] という追加の Bin があります。たまたま、Bin[0] と Bin[n/2] は常に複素数コンポーネント 0 を持つため、Bin[0] を格納するために A[0].realp が使用され、Bin[ を格納するために A[0].imagp が使用されます。 n/2]
そして、各複素数の大きさは、その周波数の周りで振動するエネルギーの量です。
したがって、ご覧のとおり、十分に細かい粒度がほとんどないため、非常に優れたピッチ検出器にはなりません。特定のビンの正確な周波数を取得するために、フレーム間の位相変化を使用して FFT ビンから正確な周波数を抽出するという狡猾なトリックがあり ます。
さて、コードに進みます：
vDSP_fft_zrip の 'ip' = 'in place' に注意してください。つまり、出力は A を上書きします ('r' は実際の入力を受け取ることを意味します)。
vDSP_fft_zrip のドキュメントを見てください。
実数データは分割複素数形式で格納され、奇数実数は分割複素数形式の虚数側に格納され、偶数実数は実数側に格納されます。
これはおそらく最も理解するのが難しいことです。プロセス全体で同じコンテナー (&A) を使用しています。そのため、最初に n 個の実数を入力します。FFT の後、n/2 の複素数を保持することになります。次に、それを逆変換に投入し、うまくいけば元の n 実数を取得します。
これで A の構造が複雑な値の設定になりました。そのため、vDSP は実数をパックする方法を標準化する必要があります。
最初に n 個の実数を生成します: 1, 2, ..., n
for (i = 0; i < n; i++)
originalReal[i] = (float) (i + 1);
次に、それらを n/2 の複雑な # として A にパックします。
// 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
これを取得するには、 A がどのように割り当てられているかを実際に確認する必要があります。ドキュメントで COMPLEX_SPLIT を調べてください。
A.realp = (float *) malloc(nOver2 * sizeof(float));
A.imagp = (float *) malloc(nOver2 * sizeof(float));
次に、事前計算を行います。
数学の素早い DSP クラス: フーリエ理論は理解するのに長い時間がかかります (私は数年前からオンとオフを繰り返してきました)
シソイドは次のとおりです。
z = exp(i.theta) = cos(theta) + i.sin(theta)
つまり、複素平面の単位円上の点です。
複素数を掛けると、角度が加算されます。したがって、z^k は単位円を飛び回り続けます。z^k は角度 k.theta で見つけることができます
z1 = 0+1i、つまり実軸から 4 分の 1 回転を選択し、z1^2 z1^3 z1^4 がそれぞれ別の 4 分の 1 回転を与え、z1^4 = 1 になることに注意してください。
z2 = -1、つまり半回転を選択します。また、z2^4 = 1 ですが、z2 はこの時点で 2 サイクルを完了しています (z2^2 も = 1)。したがって、z1 を基本周波数、z2 を第 1 高調波と考えることができます。
同様に、z3 = 「1 回転の 4 分の 3」ポイント、つまり -i は正確に 3 サイクルを完了しますが、実際には毎回 3/4 前進することは、毎回 1/4 後退することと同じです。
つまり、z3 は単に z1 ですが、反対方向です -- これはエイリアシングと呼ばれます
全波を保持するために 4 つのサンプルを選択したため、z2 は意味のある最も高い周波数です。
- z0 = 1+0i、z0^(なんでも)=1、これがDCオフセット
任意の 4 ポイント信号を z0 z1 と z2 の線形結合として表現でき ます。つまり、これらの基底ベクトルに投影しています。
しかし、「シソイドに信号を投影することはどういう意味ですか?」と尋ねるのを聞いています。
このように考えることができます: 針はシソイドの周りを回転するので、サンプル k では、針は方向 k.theta を指しており、長さは signal[k] です。シソイドの周波数と正確に一致する信号は、結果の形状をある方向に膨らませます。したがって、すべての貢献を合計すると、強力な合成ベクトルが得られます。周波数がほぼ一致する場合、膨らみは小さくなり、円の周りをゆっくりと移動します。周波数と一致しない信号の場合、寄与は互いに打ち消し合います。
http://complextoreal.com/tutorials/tutorial-4-fourier-analysis-made-easy-part-1/ は、直感的な理解を得るのに役立ちます。
しかし要点は次のとおりです。1024 個のサンプルを {z0,...,z512} に射影することを選択した場合、z0 から z512 までを事前計算することになり、それがこの事前計算ステップです。
これを実際のコードで実行している場合は、アプリのロード時にこれを 1 回実行し、アプリが終了するときに補完リリース関数を 1 回呼び出す必要があることに注意してください。何度もやらないでください。費用がかかります。
// 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);
log2n をたとえば 8 に設定すると、これらの事前計算された値を解像度 <= 2^8 を使用する任意の fft 関数に投入できることに注意してください。したがって、(究極のメモリ最適化が必要でない限り) 必要な最高解像度のセットを 1 つ作成し、それをすべてに使用してください。
次に、事前に計算したものを使用して、実際の変換を行います。
vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_FORWARD);
この時点で、A には n/2 の複素数が含まれます。最初の 1 つだけが実際には 2 つの実数 (DC オフセット、ナイキスト #) であり、複素数を装います。ドキュメントの概要では、このパッキングについて説明しています。それは非常にきちんとしています - 基本的に、変換の (複雑な) 結果を (実際の、しかし奇妙にパッケージ化された) 入力と同じメモリ フットプリントにパックすることができます。
vDSP_fft_zrip(setupReal, &A, stride, log2n, FFT_INVERSE);
そして再び... A から元の配列をアンパックする必要があります。次に、最初の配列が正確に戻っていることを確認するために比較し、事前に計算されたボビンを解放して完了です!
ちょっと待って！開梱する前に、最後に行う必要があることが 1 つあります。
// 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);
于 2010-11-04T07:41:27.783 に答える
26
実際の例を次に示します。Accelerate の vDSP fft ルーチンを使用してリモート IO オーディオ ユニットの入力で自己相関を行う C++ のスニペット。このフレームワークの使用はかなり複雑ですが、ドキュメントはそれほど悪くありません。
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);
}
}
于 2010-08-20T21:09:34.633 に答える
14
Apple の FFT フレームワークは高速だと言いますが... 正確なピッチ検出を行うには、FFT がどのように機能するかを知る必要があります (つまり、連続する各 FFT で位相差を計算して、正確なピッチを検出する必要があります)。最も優勢なビン）。
参考になるかわかりませんが、チューナー アプリ (musicianskit.com/developer.php) から Pitch Detector オブジェクトをアップロードしました。ダウンロード用の xCode 4 プロジェクトの例もあります (実装がどのように機能するかを確認できます)。
FFT 実装の例をアップロードする作業を行っています。そのため、しばらくお待ちください。更新されたら更新します。
ハッピーコーディング！
于 2012-06-14T23:05:03.860 に答える
4
別の実際の例を次に示します: https://github.com/krafter/DetectingAudioFrequency
于 2014-04-28T10:42:38.657 に答える