kopia lustrzana https://github.com/Dsplib/libdspl-2.0
added doc for psd_bartlett
rodzic
683f419847
commit
5be973b36b
119
dspl/src/psd.c
119
dspl/src/psd.c
|
@ -33,7 +33,118 @@
|
||||||
|
|
||||||
#endif
|
#endif
|
||||||
#ifdef DOXYGEN_RUSSIAN
|
#ifdef DOXYGEN_RUSSIAN
|
||||||
|
/*! ****************************************************************************
|
||||||
|
\ingroup PSD_GROUP
|
||||||
|
\fn int psd_bartlett(double* x, int n, int nfft,
|
||||||
|
fft_t* pfft, double fs,
|
||||||
|
int flag, double* ppsd, double* pfrq)
|
||||||
|
\brief Непараметрическая оценка спектральной плотности мощности (СПМ)
|
||||||
|
вещественного сигнала методом Бартлетта.
|
||||||
|
|
||||||
|
Функция рассчитывает спектральную плотность мощности \f$ X(f) \f$
|
||||||
|
выборки сигнала длительности \$n \$ отсчетов методом Бартлетта:
|
||||||
|
\f[
|
||||||
|
X(f) = \frac{1}{N F_s } \sum_{p = 0}^{P-1}\left| \sum_{m = 0}^{n_{FFT}-1}
|
||||||
|
x(m+p \cdot n_{\text{FFT}}) \exp
|
||||||
|
\left( -j 2\pi f m \right) \right|^2,
|
||||||
|
\f]
|
||||||
|
где \f$ F_s \f$ -- частота
|
||||||
|
дискретизации (Гц), \f$P = n/n_{\text{FFT}}\f$ -- количество сегментов
|
||||||
|
смещений выборки сигналов размера \f$n_{FFT}\f$.
|
||||||
|
|
||||||
|
|
||||||
|
При использовании \f$n_{FFT} = n\f$ оценка Бартлетта переходит
|
||||||
|
в стандартную периодограмму.
|
||||||
|
|
||||||
|
Расчет спектральной плотности мощности ведется при помощи алгоритмов быстрого
|
||||||
|
преобразования Фурье, для дискретной сетки частот от 0 Гц до \f$ F_s \f$ Гц
|
||||||
|
(по умолчанию), или от \f$-F_s /2 \f$ до \f$F_s /2 \f$, если установлен флаг
|
||||||
|
расчета двусторонней СПМ.
|
||||||
|
|
||||||
|
\note Метод Бартлетта возвращает асимптотически несмещенную,
|
||||||
|
состоятельную оценку СПМ (уровень флуктуаций шумовой СПМ
|
||||||
|
уменьшается с ростом длины выборки `n` при фиксированной `nfft`).
|
||||||
|
|
||||||
|
\param[in] x
|
||||||
|
Указатель на входной вектор вещественного сигнала \f$x(m)\f$,
|
||||||
|
\f$ m = 0 \ldots n-1 \f$. \n
|
||||||
|
Размер вектора `[n x 1]`. \n \n
|
||||||
|
|
||||||
|
\param[in] n
|
||||||
|
Размер вектора входного сигнала.
|
||||||
|
Также размер выходного вектора СПМ и
|
||||||
|
вектора частоты также равны `n`.\n\n
|
||||||
|
|
||||||
|
|
||||||
|
\param[in] nfft
|
||||||
|
Размер сегмента.\n
|
||||||
|
Размер выходного вектора СПМ, и соответствующего ей вектора частоты.\n\n
|
||||||
|
|
||||||
|
|
||||||
|
\param[in] pfft
|
||||||
|
Указатель на структуру \ref fft_t. \n
|
||||||
|
Указатель может быть `NULL`. В этом случае объект структуры будет
|
||||||
|
создан внутри функции и удален перед завершением.\n
|
||||||
|
Если предполагается многократный вызов функции, то рекомендуется создать
|
||||||
|
объект \ref fft_t и передавать в функцию, чтобы не
|
||||||
|
создавать его каждый раз. \n\n
|
||||||
|
|
||||||
|
\param[in] fs
|
||||||
|
частота дискретизации выборки исходного сигнала (Гц). \n\n
|
||||||
|
|
||||||
|
\param[in] flag
|
||||||
|
Комбинация битовых флагов, задающих режим расчета:
|
||||||
|
\verbatim
|
||||||
|
DSPL_FLAG_LOGMAG - СПМ считать в логарифмическом масштабе в единицах дБ/Гц
|
||||||
|
DSPL_FLAG_PSD_TWOSIDED - двусторонняя СПМ (от -Fs/2 до Fs/2)
|
||||||
|
\endverbatim
|
||||||
|
|
||||||
|
\param[in, out] ppsd
|
||||||
|
Указатель на вектор СПМ рассчитанных по входному сигналу $x$. \n
|
||||||
|
Размер вектора `[n x 1]`. \n
|
||||||
|
Память должна быть выделена. \n\n
|
||||||
|
|
||||||
|
\param[in, out] pfrq
|
||||||
|
Указатель на вектор частоты, соответствующей
|
||||||
|
значениям рассчитанного вектора СПМ. \n
|
||||||
|
Размер вектора `[n x 1]`. \n
|
||||||
|
Указатель может быть `NULL`,в этом случае вектор частоты не
|
||||||
|
рассчитывается и не возвращается. \n\n
|
||||||
|
|
||||||
|
\return
|
||||||
|
`RES_OK` если расчет произведен успешно. \n
|
||||||
|
В противном случае \ref ERROR_CODE_GROUP "код ошибки". \n \n
|
||||||
|
|
||||||
|
Пример оценок СПМ методом Бартлетта:
|
||||||
|
|
||||||
|
\include psd_bartlett_test.c
|
||||||
|
|
||||||
|
Программа производит расчет СПМ сигнала, состоящего из двух гармоник на
|
||||||
|
фоне белого гауссова шума. Расчет ведется по выборкe длины 8192 отсчета
|
||||||
|
при длине сегмента `nfft` 128, 1024 и
|
||||||
|
8192 отсчетов.
|
||||||
|
|
||||||
|
Рассчитанные СПМ выводятся на графики:
|
||||||
|
|
||||||
|
`n = 8192, nfft = 8192`:
|
||||||
|
\image html psd_bartlett_8192.png
|
||||||
|
|
||||||
|
`n = 8192, nfft = 1024`:
|
||||||
|
\image html psd_bartlett_1024.png
|
||||||
|
|
||||||
|
`n = 8192, nfft = 128`:
|
||||||
|
\image html psd_bartlett_1024.png
|
||||||
|
|
||||||
|
|
||||||
|
Можно видеть, что метод Бартлетта позволяет снизить
|
||||||
|
уровень флуктуация шума с увеличением количества сегментов.
|
||||||
|
Однако наблюдается эффект растекания спектра, который существенно ухудшает
|
||||||
|
динамический диапазон анализа.
|
||||||
|
|
||||||
|
Для более качественной оценки СПМ смотри функцию \ref psd_welch
|
||||||
|
|
||||||
|
\author Бахурин Сергей www.dsplib.org
|
||||||
|
***************************************************************************** */
|
||||||
#endif
|
#endif
|
||||||
int DSPL_API psd_bartlett(double* x, int n, int nfft,
|
int DSPL_API psd_bartlett(double* x, int n, int nfft,
|
||||||
fft_t* pfft, double fs,
|
fft_t* pfft, double fs,
|
||||||
|
@ -815,7 +926,7 @@ DSPL_FLAG_PSD_TWOSIDED - двусторонняя СПМ (от -Fs/2 до Fs/2)
|
||||||
|
|
||||||
Пример периодограммных оценок СПМ для различной длины выборки сигнала:
|
Пример периодограммных оценок СПМ для различной длины выборки сигнала:
|
||||||
|
|
||||||
\include psd_welch_cmplx_test.c
|
\include psd_welch_test.c
|
||||||
|
|
||||||
Программа производит расчет СПМ сигнала, состоящего из двух комплексных
|
Программа производит расчет СПМ сигнала, состоящего из двух комплексных
|
||||||
гармоник на фоне белого гауссова шума.
|
гармоник на фоне белого гауссова шума.
|
||||||
|
@ -824,13 +935,13 @@ DSPL_FLAG_PSD_TWOSIDED - двусторонняя СПМ (от -Fs/2 до Fs/2)
|
||||||
Рассчитанные СПМ выводятся на графики:
|
Рассчитанные СПМ выводятся на графики:
|
||||||
|
|
||||||
`nfft = 8192, noverlap = 4096`:
|
`nfft = 8192, noverlap = 4096`:
|
||||||
\image html psd_welch_cmplx_8192.png
|
\image html psd_welch_8192.png
|
||||||
|
|
||||||
`nfft = 1024, noverlap = 512`:
|
`nfft = 1024, noverlap = 512`:
|
||||||
\image html psd_welch_cmplx_1024.png
|
\image html psd_welch_1024.png
|
||||||
|
|
||||||
`nfft = 256, noverlap = 128`:
|
`nfft = 256, noverlap = 128`:
|
||||||
\image html psd_welch_cmplx_256.png
|
\image html psd_welch_256.png
|
||||||
|
|
||||||
Можно видеть, что уменьшение `nfft` при фиксированной длительности сигнала
|
Можно видеть, что уменьшение `nfft` при фиксированной длительности сигнала
|
||||||
позволяет уменьшить флуктуации шума и делает оценку состоятельной.
|
позволяет уменьшить флуктуации шума и делает оценку состоятельной.
|
||||||
|
|
Ładowanie…
Reference in New Issue