Python — Расчет PSD

Я пытаюсь создать спектрограмму, как в matplotlib и в этой библиотеке они возвращают PSD (Power Spectra Density), а не абсолютную величину и т. д.

Я хочу внедрить PSD в свой собственный код (написанный на C ++). Я искал исходный код для matplotlib это выглядит примерно так:

  result, windowVals = apply_window(result, window, axis=0,
return_window=True)
result = np.fft.fft(result, n=pad_to, axis=0)[:numFreqs, :]

result = np.conjugate(result) * result

# calculations for PSD

result /= (np.abs(windowVals)**2).sum()

result[1:-1] *= scaling_factor

Ссылка может быть найдена: Вот

Теперь я хотел бы сделать то же самое, я дошел до той части, где я применил окна и выполнил DFT. Оттуда, где я думаю, что я иду не так.

Вот мой код для расчета PSD.

std :: vector CalculatePSD (std :: vector &Vals)
{

std::vector<double> hanning = getHanningWindow(vals.size());

std::vector<double> vals2(vals.size());

for(unsigned i=0; (i < vals.size()); i++)
{
double v = vals[i].re * (-1 * vals[i].im);
vals2[i] = v * v;
}

double vi = 0.0;
for(unsigned i=0; (i < vals2.size()); i++)
{
vi += abs(hanning[i]);
}

for(unsigned i=0; (i < vals2.size()); i++)
{
vals2[i] /= vi;

}

for(unsigned i=0; (i < vals2.size()); i++)
{
vals2[i] *= 2;
}
return vals2;
}

Мои результаты, которые я получаю:

[  2.66926000e-10   4.71270000e-10   3.25024000e-10   1.86008000e-10
8.68276000e-11   2.53098000e-11   7.92737000e-13   5.29274000e-12
1.90057000e-11   3.27736000e-11   4.42093000e-11   4.60725000e-11
4.98549000e-11   6.08991000e-11   5.59033000e-11   2.26625000e-11
7.17159000e-13   1.65150000e-12   1.31530000e-12   7.86863000e-14
1.59211000e-13   2.80960000e-13   1.25471000e-13   1.35225000e-11
1.09812000e-10   3.49411000e-10   5.96210000e-10   6.43539000e-10
4.85667000e-10   2.46446000e-10   7.61480000e-11   7.70735000e-12
3.79924000e-13   4.16633000e-13   9.22683000e-12   6.23747000e-11
1.80818000e-11   5.13819000e-11   7.20294000e-10   2.46505000e-09
4.42844000e-09   4.77968000e-09   3.06048000e-09   1.12040000e-09
1.15313000e-10   4.78101000e-11   2.96048000e-10   3.93216000e-10
2.66339000e-10   1.01348000e-10   3.13012000e-11   3.59365000e-11
3.00758000e-11   2.44970000e-11   4.19013000e-12   1.92963000e-12
7.62176000e-13   9.83650000e-14   4.47749000e-15   3.29484000e-19
7.46089000e-14   3.99184000e-12   4.75458000e-11   3.36518000e-10
9.34579000e-10   9.14683000e-10   1.34570000e-10   6.47235000e-11
1.21893000e-10   1.44874000e-11   5.37149000e-15   1.93477000e-11
1.15593000e-10   1.27684000e-10   2.10146000e-11   9.34294000e-13
8.85961000e-12   6.35360000e-12   3.88437000e-13   2.75095000e-13
1.07261000e-11   9.35213000e-12   4.16280000e-11   1.52995000e-11
9.86730000e-12   2.04171000e-11   9.69115000e-12   1.81067000e-12
1.64522000e-13   4.01619000e-12   2.81457000e-13   5.86617000e-13
1.33937000e-11   3.16438000e-11   8.83417000e-11   1.24155000e-10
9.36973000e-11   1.89000000e-12   4.76718000e-14   3.46776000e-13
3.88868000e-12   1.01426000e-12   7.87006000e-13   1.02698000e-10
6.80669000e-12   3.22014000e-12   9.92309000e-12   1.17853000e-10
1.20818000e-11   2.20125000e-12   8.78943000e-12   3.34323000e-10
4.28139000e-10   3.57678000e-10   1.57808000e-10   2.05267000e-11
8.98399000e-11   1.12894000e-11   3.21320000e-14   7.77191000e-12
5.91681000e-13   9.92243000e-16   2.10894000e-11   2.19397000e-12
1.14148000e-12   3.41732000e-12   1.81439000e-12   1.46689000e-11]

Принимая во внимание, что результаты в Python:

[  6.44554713e-04   2.26979569e-02   1.48395306e-02   1.39560086e-02
1.70585613e-02   4.24042116e-04   4.10722082e-04   1.77314474e-02
5.48046037e-03   6.86724979e-03   1.33342952e-02   5.45918807e-04
1.42011959e-06   8.15283041e-03   3.02976247e-02   2.95310636e-02
2.69222586e-02   2.70161073e-04   4.27988811e-04   8.22069685e-03
1.14550280e-03   5.94684341e-03   5.03412155e-03   2.39065158e-04
1.88851349e-03   1.63618611e-02   1.02155767e-02   5.56409334e-03
2.03783039e-02   1.30646965e-03   7.83925381e-03   6.58153969e-04
8.58222471e-05   4.90329132e-03   9.27321780e-03   2.18878971e-02
7.80419597e-03   1.65506496e-05   2.12233732e-03   3.48564618e-02
3.04324943e-02   1.14097124e-02   1.83163044e-02   5.53528648e-04
1.72024876e-03   1.05496508e-02   1.22350425e-02   6.81764861e-03
2.18181750e-02   1.25305967e-04   3.45533908e-04   2.52806605e-02
2.79032703e-02   3.30741745e-02   8.92045889e-03   1.43861624e-04
1.37729407e-03   4.40048633e-02   4.43466583e-02   3.21348174e-02
1.97845126e-02   4.76052263e-05   1.90059116e-03   1.36124930e-02
3.08483724e-02   3.18817777e-02   4.22224299e-02   6.48991341e-04
3.37298579e-04   1.92796091e-02   3.26384995e-02   7.44582047e-03
2.75911372e-02   8.24143406e-05   6.19298800e-04   2.18909904e-02
7.85534253e-03   1.35622242e-02   1.64534364e-02   4.24610034e-07
1.25296965e-04   3.85049720e-03   1.56208315e-02   1.51067447e-02
1.15560295e-02   1.91845524e-02   1.51484986e-02   3.68090803e-05
4.56878093e-04   1.32583767e-02   2.67413477e-02   2.12116190e-02
1.40731136e-02   7.46782595e-06   7.56130481e-04   1.11894743e-02
3.82556474e-02   2.20488800e-02   1.18449472e-02   6.41610843e-05
8.93214478e-04   1.44705708e-02   8.95544599e-03   8.24627650e-03
1.54125088e-02   3.82922435e-07   1.21567170e-03   3.66207393e-02
2.52421164e-02   2.79258696e-02   2.42711875e-02   4.41070028e-04
9.18506931e-04   2.29391748e-02   2.93676503e-06   3.51546485e-02
1.53622376e-03   3.93588210e-05   3.35935113e-04   1.74232319e-02
1.89744096e-02   1.02178421e-02   1.54763125e-02   7.24992746e-05
3.46909205e-04   2.41130633e-02   1.96140922e-02   1.94479820e-02
1.02416629e-02   1.08582494e-04   5.91329398e-04   8.53889890e-03
8.05471223e-03   1.78734494e-02   2.12358089e-02   2.74402258e-02
1.86277585e-02   2.97114647e-06   1.58242458e-03   1.77131478e-02
1.03301989e-03   1.17236867e-02   2.70723000e-02   2.45157582e-04
1.70524416e-04   1.12361676e-03]

Я, вероятно, что-то не так с преобразованием кода, так как я плохо знаю Python. Если кто-то определит, где я могу пойти не так, это очень мне поможет.

РЕДАКТИРОВАТЬ:

    complex<double> foo[vals.size()];
for(unsigned i=0; (i < vals.size()); i++)
{
foo[i] = complex<double>(vals[i].re, vals[i].im);
}

std::vector<double> vals2(vals.size());

for(unsigned i=0; (i < vals2.size()); i++)
{
double d =(std::conj(foo[i])*std::real(foo[i]));
}

0

Решение

Ваша функция возвращает vals (изменяя его по ссылке), пока вы, кажется, сохраняете результат в vals2, которая является локальной переменной, которая удаляется в конце функции.

РЕДАКТИРОВАТЬ

Да, похоже, что ошибка может быть в сопряжении — почему бы не использовать функциональность std :: complex, воспроизводящую код Python, а не изобретать колесо:

for(unsigned i=0; (i < vals.size()); i++)
{
vals2[i]= (vals[i].conj()*vals[i]).re();
}
0

Другие решения

Других решений пока нет …

По вопросам рекламы [email protected]