bool BilateralFilter(const cv::Mat& color, cv::Mat& depth, double sigma_w, double sigma_c, int w) {
    if (color.channels() != 4) {
        return false;
    }

    // Ensure depth is of type CV_32F
    if (depth.type() != CV_32F) {
        return false;
    }

    cv::Mat colorDouble;
    color.convertTo(colorDouble, CV_32FC4); // Convert to double (RGBA)
    cv::Mat depthDouble;
    depth.convertTo(depthDouble, CV_32F); // Convert to double

    int inputHeight = color.rows;
    int inputWidth = color.cols;
    
    cv::Mat result = cv::Mat::zeros(inputHeight, inputWidth, CV_32F);

    // Precompute the spatial Gaussian kernel
    cv::Mat spatialKernel = cv::Mat(2 * w + 1, 2 * w + 1, CV_32F);
    for (int mx = -w; mx <= w; ++mx) {
        for (int my = -w; my <= w; ++my) {
            spatialKernel.at<float>(mx + w, my + w) = std::exp(-(mx * mx + my * my) / (2 * sigma_w * sigma_w));
        }
    }

    for (int i = 0; i < inputHeight; ++i) {
        for (int j = 0; j < inputWidth; ++j) {
            // 提取局部区域
            int iMin = std::clamp(i - w, 0, inputHeight - 1);
            int iMax = std::clamp(i + w, 0, inputHeight - 1);
            int jMin = std::clamp(j - w, 0, inputWidth - 1);
            int jMax = std::clamp(j + w, 0, inputWidth - 1);

            cv::Rect roi(jMin, iMin, jMax - jMin + 1, iMax - iMin + 1);
            cv::Mat colorSec = colorDouble(roi);
            cv::Mat depthSec = depthDouble(roi);
            if (cv::countNonZero(depthSec) == 0) {
                result.at<float>(i, j) = 0;
                continue;
            }
           
            // Split the RGBA channels
            std::vector<cv::Mat> channels(4);
            cv::split(colorSec, channels);

            // Compute Gaussian range weights using RGB channels (ignoring the alpha channel)
            cv::Mat dR = channels[0] - colorDouble.at<cv::Vec4f>(i, j)[0];
            cv::Mat dG = channels[1] - colorDouble.at<cv::Vec4f>(i, j)[1];
            cv::Mat dB = channels[2] - colorDouble.at<cv::Vec4f>(i, j)[2];
            
            // Calculate the squared difference
            cv::Mat rangeWeight = -(dR.mul(dR) + dG.mul(dG) + dB.mul(dB)) / (2 * sigma_c * sigma_c);
            
            // Apply the exponential function element-wise
            cv::Mat range = cv::Mat(rangeWeight.size(), rangeWeight.type());
            for (int y = 0; y < rangeWeight.rows; ++y) {
                for (int x = 0; x < rangeWeight.cols; ++x) {
                    range.at<float>(y, x) = std::exp(rangeWeight.at<float>(y, x));
                }
            }

            // Calculate bilateral filter response
            cv::Mat spatialSubsection = spatialKernel(cv::Rect(jMin - j + w, iMin - i + w, roi.width, roi.height));
            cv::Mat depthWeight = range.mul(spatialSubsection);
            // Apply the condition depthSec > 0 to depthWeight
            for (int y = 0; y < depthSec.rows; ++y) {
                for (int x = 0; x < depthSec.cols; ++x) {
                    if (depthSec.at<float>(y, x) <= 0) {
                        depthWeight.at<float>(y, x) = 0;
                    }
                }
            }

            cv::Mat depthSum = depthSec.mul(depthWeight);
            
            float sumDepthWeight = cv::sum(depthWeight)[0];
            result.at<float>(i, j) = cv::sum(depthSum)[0] / sumDepthWeight;
        }
    }

    // Convert result back to CV_32F and store in depth
    result.convertTo(depth, CV_32F);
    
    return true;
}  理解这段代码，并寻找优化该代码的方案，最终给出优化代码