0

我正在学习 FBP,在过去的 1 个月里,我一直坚持使用反投影算法。到目前为止,我已经完成了以下步骤。

  • 创建正弦图(Radon 变换)

  • 创建过滤器(汉宁过滤器:我使用汉宁窗口函数来实现过滤器)我知道 Ram 过滤器效果最好,但我想研究不同的高通过滤器将如何工作。

  • 快速傅里叶变换。(使用dft函数)

  • 应用过滤器

  • 傅里叶逆变换(使用dft函数)

  • 将滤波后的正弦图通过反投影算法。我在这里分享两个功能:

  • 传统的反投影算法(来自LOC:29

Mat backprojectionTraditional(Mat inversefft) {
    
    Mat reconstruction(inversefft.size(), CV_32F);
    
    int numOfAngles = 180;
    int dtheta = 180 / numOfAngles;
    for (int x = 0; x < reconstruction.size().height; x++) {
        double delta_t = M_PI / 180;
        for (int y = 0; y < reconstruction.size().width; y++) {
            reconstruction.at<float>(x, y) = 0.0;
            for (int theta = 0; theta < inversefft.size().width; theta += dtheta) {
                int s = (x - 0.5 * inversefft.size().height) * sin(theta * delta_t) +
                    (y - 0.5 * inversefft.size().height) * cos(theta * delta_t) + 0.5 * inversefft.size().height;
                if (s > 0 && s < inversefft.size().height) {
                    reconstruction.at<float>(x, y) += inversefft.at<float>(s, theta);
                }
            }
            if (reconstruction.at<float>(x, y) < 0)reconstruction.at<float>(x, y) = 0;
        }
    }
    rotate(reconstruction, reconstruction, ROTATE_90_CLOCKWISE);
    return reconstruction;
}
Mat backprojectionIterative(Mat inversefft) 
{
    Mat iradon(inversefft.size(), CV_32F);
    Mat projection(inversefft.size(), CV_32F);
    float colsum;
    for (int z = 0; z < 180; z++) {
        projection = imrotate(inversefft, z);
        //imshow("Projected image", projection);
        //waitKey(0);
        
        for (int i = 0; i < inversefft.size().width; i++) {
            colsum = 0;
            for (int j = 0; j < inversefft.size().height; j++) {
                colsum += projection.at<float>(j, i);
                //cout << colsum << endl;

            }
            for (int k = 0; k < inversefft.size().height; k++) {
                //cout << "Break" << endl;
                //cout << colsum << endl;
                //cout << iradon.at<float>(k, i) << endl;
                iradon.at<float>(k, i) += colsum;
                //cout << iradon.at<float>(k, i) << endl;
            }
        }
        
        //iradon = imrotate(iradon, 90);
        //normalize(iradon, iradon, 0, 1, NORM_MINMAX, CV_32F);
        //imshow("Reconstructed image", iradon);
        //waitKey(0);
    }
    return iradon;
}

谁能告诉我如何实现反投影算法以及我的错误在哪里?

这是测试数据

#PS - 这是我的第一个问题,如果我不遵守社区标准,我会尽快解决。谢谢你。

编辑 - 由于大多数观众不理解问题陈述,我正在上传这两个函数的结果。
应用backprojectionTraditional功能
后 应用backprojectionIterative功能后

应用反投影后,输出应与输入相同或接近。

4

0 回答 0