Почему это лучше для использования'! =", чем' <' в векторном цикле? (C++)

1) Как вы правильно поняли, требование к масштабированию связано с дискретным преобразованием Фурье. Лучший способ получить это - вычислить деконволюцию двух унифицированных единичных сигналов. Их DFT равен n 0 0 0 ...., где n - количество точек DFT. Следовательно, отношение r_f равно 1 0 0 0 0, и его обратное БПФ, вычисленное с помощью np.fft.ifft(), равно 1 / n 1 / n 1 / n ... Правильный сигнал, полученный в результате деконволюции, должен был быть 1 / T 1 / T 1 / Т ..., где Т = 10. это длина кадра.

Как следствие, правильное масштабирование для выполнения деконволюции составляет n/T= len(r_f)/10.

r_if=r_if*len(r_if)/10.

2) Деконволюционный сигнал транслируется на половину периода. Это связано с тем, что ядро ​​Гаусса центрируется по центру кадра. Просто сдвиньте ядро ​​на половину периода, и проблема решена. Для этого можно применить функцию np.fft.fftshift():

f_yb = np.fft.fft(np.fft.fftshift(y_b))    #fft the filter

РЕДАКТИРОВАТЬ: Чтобы исследовать причину перевода, давайте сосредоточимся на случае, когда ядро ​​деконволюции является очень узким гауссовым распределением, почти соответствующим Распределение Дирака. Ваш входной сигнал представляет собой гауссову кривую с центром в нуле, а кадр дискретизируется между -5 и 5. Аналогично, ядро ​​деконволюции представляет собой дираковский центр с центром в нуле. Как следствие, деконволюция сигнала должна быть идентична входному сигналу: гауссова кривая с центром в нуле. Тем не менее, DFT, реализованный в FFTW и, следовательно, np.fft.fft() , вычисляется как кадр, начинающийся с 0 и заканчивающийся на 10 , который выбирается в точках 10j / n, где j находится в [0..n- 1], частоты в пространстве Фурье равны k / 10, где k в [0..n / 2, -n / 2 + 1 ..- 1]. Как следствие, этот ДПФ воспринимает ваш сигнал как гауссиан с центром в 5, а ядро ​​деконволюции - как дираковский с центром в 5. Свертка функции f (t) с дельтой Дирака (t-t0) с центром в точке t0 является просто переведенная функция f (t-t0). Следовательно, результатом деконволюции, вычисленной в np.fft.fft(), является входной сигнал, переведенный на половину периода. Поскольку входной сигнал центрируется на 0 в кадре [-5,5], выходной сигнал, вычисленный с помощью np.fft.fft(), центрируется на -5 (или эквивалентно 5 из-за периодичности). Сдвиг ядра устраняет несоответствие между тем, что мы думаем о кадре как [-5 5], и np.fft.ifft() обрабатываем его так, как если бы это было [0 10].

Фильтры часто предназначены для уменьшения влияния высокочастотных шумов. Таким образом, деконволюция вызывает потенциальное увеличение высокочастотного шума. Скрининг частот, как вы сделали, является потенциальным решением. Обратите внимание, что это точно эквивалентно свертке сигнала с определенным фильтром!

В диапазоне томографической реконструкции алгоритм отфильтрованного обратного проецирования требует применения линейного фильтра, который значительно раздувает высокочастотный шум. Здесь предложен фильтр Винера : этот тип фильтра может быть спроектирован так, чтобы минимизировать среднеквадратичную ошибку на деконволюционном сигнале, учитывая SNR извилистого сигнала. Тем не менее, требуется некоторое предположение относительно спектральной плотности мощности сигнала и шума.