有限数のサンプルで均一にサンプリングされた離散信号に対してSincダウンサンプリング（DFTダウンサンプリング）を行う適切な方法
与えられた信号 $ \left\{ x [ 0 ], x [ 1 ], ..., x [ N - 1 ] \right\} $ 周波数領域でダウンサンプリングする正しい方法は何ですか（Sinc補間）？
回答
周波数内挿（DFTドメイン）
実装はよく知られています。MATLABでは、次のようになります。
if(numSamplesO > numSamples)
% Upsample
halfNSamples = numSamples / 2;
if(mod(numSamples, 2) ~= 0) % Odd number of samples
vXDftInt = interpFactor * [vXDft(1:ceil(halfNSamples)); zeros(numSamplesO - numSamples, 1, 'like', vXDft); vXDft((ceil(halfNSamples) + 1):numSamples)];
else % Even number of samples -> Special Case
vXDftInt = interpFactor * [vXDft(1:halfNSamples); vXDft(halfNSamples + 1) / 2; zeros(numSamplesO - numSamples - 1, 1, 'like', vXDft); vXDft(halfNSamples + 1) / 2; vXDft((halfNSamples + 2):numSamples)];
end
else
% Downsample
halfNSamples = numSamplesO / 2;
if(mod(numSamples, 2) ~= 0) % Odd number of samples
vXDftInt = interpFactor * [vXDft(1:ceil(halfNSamples)); vXDft((numSamples - floor(halfNSamples) + 1):numSamples)];
else % Even number of samples -> Special Case
vXDftInt = interpFactor * [vXDft(1:halfNSamples); vXDft(halfNSamples + 1) / 2; vXDft((numSamples - halfNSamples + 2):numSamples)];
end
end
したがって、ここでは2つのケースを処理します。
- アップサンプル出力のサンプル数に一致するよう
に、DFTの中央部分にゼロサンプルを追加します（numSamplesO）。
入力サンプル数（numSamples）が偶数の場合に対応します。その場合、ナイキストサンプルを分割します（$ X \left[ N / 2 \right] $）2に $ N $ はサンプルの入力数です。 - ダウン
サンプル出力のサンプル数に一致するように、DFTの中央部分のサンプルを削除します（numSamplesO）。
サンプルの出力数（numSamplesO）が偶数の場合に対応します。その場合、を分割してナイキストサンプルにします（$ X \left[ M / 2 \right] $）2に $ M $ はサンプルの出力数です。
問題は、なぜこのようにするのかということです。なぜ補間係数interpFactor？の分割係数はどこにありますか$ 0.5 $から来る？
DFTを覚えておく必要があると答えるには、基本的に離散フーリエ級数（DFS）です。
これは、最も重要な仮定は、データが時間領域と周波数領域の両方で周期的であることを意味します。
さて、DFTは基本的にDFSであるため、その周期内で信号を補間する自然な方法は、フーリエ級数を使用することです。
詳細に入る前に、インデックスの値を定義するために使用される2セットの整数を定義しましょう。
$$ \begin{aligned} \mathcal{K}_{DFS}^{N} & = \left\{- \left\lceil \frac{N - 1}{2} \right\rceil, - \left\lceil \frac{N - 1}{2} \right\rceil + 1, \ldots, -1, 0, 1, \ldots, \left\lceil \frac{N - 1}{2} \right\rceil - 1, \left\lceil \frac{N - 1}{2} \right\rceil \right\} \\ \mathcal{K}_{DFT}^{N} & = \left\{- \left\lceil \frac{N - 1}{2} \right\rceil, - \left\lceil \frac{N - 1}{2} \right\rceil + 1, \ldots, -1, 0, 1, \ldots, \left\lceil \frac{N - 1}{2} \right\rceil - 1, \left\lfloor \frac{N - 1}{2} \right\rfloor \right\} \\ \end{aligned} $$
これは、最大帯域幅が $ \frac{1}{2 T} $ のサンプリング定理によってサンプリングされた $ t \in \left[ 0, N T \right) $ どこ $ T $ はサンプリング期間であり、 $ P = N T $ 機能期間です：
$$ \begin{aligned} x \left( t \right) {\Big|}_{t = n T} & = \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lceil \frac{N - 1}{2} \right\rceil} {c}_{k} {e}^{ j 2 \pi \frac{k t}{P} } && \text{By Fourier Series} \\ & = \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lceil \frac{N - 1}{2} \right\rceil} {c}_{k} {e}^{ j 2 \pi \frac{k t}{N T} } && \text{By the period of the function / series} \\ & = \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lceil \frac{N - 1}{2} \right\rceil} {c}_{k} {e}^{ j 2 \pi \frac{k n}{N} } && \text{Setting $ t = n T $} \\ & = \frac{1}{N} \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lfloor \frac{N - 1}{2} \right\rfloor} X \left[ k \right] {e}^{ j 2 \pi \frac{k n}{N} } && \text{The DFT} \end{aligned} $$
上記の式は偶数の場合でも機能します $ N = 2 l, \; l \in \mathbb{N} $ そして奇妙な場合のために $ N = 2 l + 1, \; l \in \mathbb{N} $。上記は、DFT係数とフーリエ級数係数の間の関係を定義しています。
$$ {c}_{k} = \begin{cases} \frac{ X \left[ k \right ] }{2 N} & \text{ if } k = \frac{N}{2} \\ \frac{ X \left[ k \right ] }{2 N} & \text{ if } k = -\frac{N}{2} \\ \frac{ X \left[ k \right ] }{N} & \text{ if } k \notin \left\{\frac{N}{2}, -\frac{N}{2} \right\} \end{cases}, \; k \in \mathcal{K}_{DFS}^{N} $$
しかし、他のサンプリングポイントを任意のセットに使用することを妨げるものは何もありません $ { \left\{ {t}_{m} \right\}}_{m = 0}^{M - 1} $ どこ $ \forall m, {t}_{m} \in \left[ 0, N T \right) $。これは$ x \left( t \right) = \frac{1}{N} \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lfloor \frac{N - 1}{2} \right\rfloor} X \left[ k \right] {e}^{ j 2 \pi \frac{k t}{N T} } $ ために $ t \in \left[ 0, N T \right) $。これは、複雑で実際の信号に対して機能します。
実際の信号の場合、$ x \left( t \right) \in \mathbb{R} $DFTのコサイン形式を使用することもできます。
$$ \begin{aligned} x \left( t \right) & = \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lceil \frac{N - 1}{2} \right\rceil} {c}_{k} {e}^{ j 2 \pi \frac{k t}{N T} } && \text{From the above} \\ & = \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lceil \frac{N - 1}{2} \right\rceil} \left| {c}_{k} \right| \cos \left( 2 \pi \frac{k t}{N T} + \angle {c}_{k} \right) && \text{Fourier series in its Cosine form} \\ & = \sum_{k = - \left\lceil \frac{N - 1}{2} \right\rceil}^{\left\lfloor \frac{N - 1}{2} \right\rfloor} \frac{\left| X \left[ k \right] \right|}{N} \cos \left( 2 \pi \frac{k t}{N T} + \angle X \left[ k \right] \right) && \text{Fourier series in its Cosine form} \\ & = \sum_{k = 0}^{\left\lfloor \frac{N - 1}{2} \right\rfloor} {\alpha}_{k} \frac{\left| X \left[ k \right] \right|}{N} \cos \left( 2 \pi \frac{k t}{N T} + \angle X \left[ k \right] \right) && \text{Using the DFT conjugate symmetry of a real signal} \end{aligned} $$
どこ $ {\alpha}_{k} = \begin{cases} 1 & \text{ if } k \in \left\{ 0, \frac{N}{2} \right\} \\ 2 & \text{ else } \end{cases} $。
したがって、ここで見たものと、それが上記のアルゴリズムにどのように関連しているかを考える必要があります。
まず、ここでの主なトリックは、インデックスが作成されたときにDFTのネイティブ形式が必要であることに注意する必要があります。$ k \in \mathcal{K}_{DFT}^{N} $。そうすれば、DFTの離散フーリエ級数（DFS）の原点への接続を簡単に確認できます。
備考：実際には、DFTは次のように定義（および計算）されます。$ k \in \left\{ 0, 1, \ldots, N - 1 \right\} $。
出力均一時間グリッドのセットを選択した場合 $ { \left\{ {t}_{m} \right\}}_{m = 0}^{M - 1} $ 形になります $ {t}_{m} = m {T}_{s} $ ここで、アップサンプリングレート（ダウンサンプリングは後で処理します） $ q = \frac{M}{N} \geq 1 $次に、グリッドを回復するためにIDFTを調べることによって何をする必要があるかは明らかです。
$$ x \left[ m \right] = \frac{1}{M} \sum_{k = 0}^{M - 1} \tilde{X} \left[ k \right] {e}^{j 2 \pi \frac{k m}{M}} = \frac{1}{M} \sum_{k = - \left\lceil \frac{M - 1}{2} \right\rceil}^{\left\lfloor \frac{M - 1}{2} \right\rfloor} \tilde{X} \left[ k \right] {e}^{j 2 \pi \frac{k m}{M}} $$
次に、これを上からの補間式と一致させる必要があります。これは線形変換であるため、$ q $定数の世話をします。また、$ \forall m, \frac{m}{M} = \frac{{t}_{m}}{N T} $ したがって、次のように設定します。
$$ \tilde{X} \left[ k \right] = \begin{cases} X \left[ k \right] & \text{ if } k \in \mathcal{K}_{DFT}^{N} \setminus \left\{ k \mid k = \frac{N}{2} \right\} \\ \frac{X \left[ k \right]}{2} & \text{ if } k = \frac{N}{2} \\ 0 & \text{ if } k \notin \mathcal{K}_{DFT}^{N} \end{cases} $$
から $ N $ DFTの周期性は、補間係数が次の均一な時間グリッドの最終補間を記述できます。 $ q $：
$$ x \left[ m \right] = \frac{q}{M} \sum_{k = 0}^{M - 1} \hat{X} \left[ k \right] {e}^{j 2 \pi \frac{k m}{M}} $$
どこ $ \hat{X} \left[ k \right] $ と定義されている：
$$ \hat{X} \left[ k \right] = \begin{cases} X \left[ k \right] & \text{ if } k \in \left\{ 0, 1, \ldots, N - 1 \right\} \setminus \left\{ \frac{N}{2} \right\} \\ \frac{X \left[ k \right]}{2} & \text{ if } k = \frac{N}{2} \\ 0 & \text{ if } k \in \left\{ N, N + 1, \ldots, M - 1 \right\} \end{cases} $$
上記のアップサンプルコードで行ったことはまさにこれです。
ダウンサンプルはどうですか？さて、コードが示すように、DFTドメインで同じ直感を使用できます。これは基本的に、フーリエ級数係数を使用した補間が、周波数領域でのディリクレ核による乗算に他ならないためです。$ \operatorname{sinc} \left( \cdot \right) $関数。これは、$ \frac{1}{2} $ジャンプの不連続性がある周波数領域で値1の長方形を乗算するため、因数分解します。実際、フーリエ級数は不連続点でのジャンプの平均値に収束します。から行くので$ 1 $ に $ 0 $、ジャンプ時の値が $ 0.5 $。
したがって、上記のダウンサンプリングおよびアップサンプリングコードは、アップサンプルの場合は入力のサンプリング頻度に応じて、ダウンサンプルの場合は出力に応じて、ディリクレ核をデータに適用します。
ダウンサンプリングする別の方法は、サンプルの出力数の整数係数にアップサンプリングすることです。次に、デシメーション（すべての...サンプルを取得）を使用してサンプルを取得します。2は、データに低レートとサンプリングレートの間の周波数にエネルギーがない場合に一致します。もしそうなら、それらは一致しません。
MATLABコードを追加します...
備考：この回答はアップサンプリングも対象としています。アップサンプリングに関する別の質問を開くことを検討するか、これを広げてください。