Hanjie's Blog

cv::Mat ColorMat(const cv::Mat &mat_float, bool draw_colorbar = false, const bool is_white_background = false, double min_val = 1, double max_val = 0, const cv::Mat &user_color = cv::Mat(), int colorbar_width = 50, int colorbar_gap = 5) {
    if (min_val > max_val) {
        cv::minMaxLoc(mat_float, &min_val, &max_val);
    }

    cv::Mat mat;
    mat_float.convertTo(mat, CV_8UC1, 255 / (max_val - min_val), -255 * min_val / (max_val - min_val));

    cv::Mat mat_show;

    if (user_color.empty()) {
        cv::applyColorMap(mat, mat_show, cv::COLORMAP_JET);
    } else {
        cv::applyColorMap(mat, mat_show, user_color);
    }

    if (is_white_background) {
        cv::Mat mask;
        cv::threshold(mat, mask, 0, 255, cv::THRESH_BINARY_INV);
        cv::Mat img_white(mat.size(), CV_8UC3, cv::Scalar(255, 255, 255));
        img_white.copyTo(mat_show, mask);
    }

    if (draw_colorbar) {
        cv::Mat color_bar_value(cv::Size(colorbar_width, mat_show.rows), CV_8UC1);
        cv::Mat color_bar;

        for (int i = 0; i < mat_show.rows; i++) {
            uchar value = 255 - 255 * float(i) / float(mat_show.rows);
            for (int j = 0; j < colorbar_width; j++) {
                color_bar_value.at<uchar>(i, j) = value;
            }
        }

        if (user_color.empty()) {
            cv::applyColorMap(color_bar_value, color_bar, cv::COLORMAP_JET);
        } else {
            cv::applyColorMap(color_bar_value, color_bar, user_color);
        }

        cv::Mat mat_colorbar_show(cv::Size(mat_show.cols + colorbar_width + colorbar_gap, mat_show.rows), CV_8UC3, cv::Scalar(255, 255, 255));

        mat_show.copyTo(mat_colorbar_show(cv::Rect(0, 0, mat_show.cols, mat_show.rows)));
        color_bar.copyTo(mat_colorbar_show(cv::Rect(mat_show.cols + colorbar_gap, 0, color_bar.cols, color_bar.rows)));

        cv::putText(mat_colorbar_show, ToStr(max_val), cv::Point(mat_show.cols + colorbar_gap, 20), cv::FONT_HERSHEY_SIMPLEX, 0.5, cv::Scalar(255, 255, 255), 1);
        cv::putText(mat_colorbar_show, ToStr(min_val), cv::Point(mat_show.cols + colorbar_gap, mat_show.rows - 10), cv::FONT_HERSHEY_SIMPLEX, 0.5, cv::Scalar(255, 255, 255), 1);

        mat_show = mat_colorbar_show;
    }

    return mat_show;
};

cv::Mat ShowColorMat(const std::string &name, const cv::Mat &mat_float, bool draw_colorbar = false, const float scale = 1, const bool is_white_background = false, double min_val = 1, double max_val = 0, const cv::Mat &user_color = cv::Mat(), int colorbar_width = 50, int colorbar_gap = 5) {
    cv::Mat mat_resize;
    cv::resize(mat_float, mat_resize, cv::Size(), scale, scale, cv::INTER_NEAREST);
    cv::Mat img = ColorMat(mat_resize, draw_colorbar, is_white_background, min_val, max_val, user_color, colorbar_width, colorbar_gap);
    cv::imshow(name, img);
    return img;
}

cv::Mat FilterMask(const float radius, const cv::Size &mask_size) {
    cv::Mat mask = cv::Mat::ones(mask_size, CV_8UC1);
    cv::circle(mask, cv::Point(mask_size.width / 2, mask_size.height / 2), radius, cv::Scalar(0), -1);
    
    // https://datahacker.rs/opencv-discrete-fourier-transform-part2/
    cv::Mat mask_float, filter_mask;
    mask.convertTo(mask_float, CV_32F);
    std::vector<cv::Mat> mask_merge = {mask_float, mask_float};
    cv::merge(mask_merge, filter_mask);
    return filter_mask;
}

bool BlurDetection(const cv::Mat &img,
                   const cv::Mat &filter_mask,
                   float &blur_value,
                   const float thresh = 10,
                   const bool debug_show = false) {
    // https://pyimagesearch.com/2020/06/15/opencv-fast-fourier-transform-fft-for-blur-detection-in-images-and-video-streams/
    // https://github.com/Qengineering/Blur-detection-with-FFT-in-C/blob/master/main.cpp
    // https://docs.opencv.org/3.4/d8/d01/tutorial_discrete_fourier_transform.html

    cv::Mat img_scale, img_gray, img_fft;
    cv::resize(img, img_scale, filter_mask.size());

    if (img_scale.channels() == 3) {
        cv::cvtColor(img_scale, img_gray, cv::COLOR_BGR2GRAY);
    } else {
        img_gray = img_scale;
    }
    img_gray.convertTo(img_fft, CV_32F);

    // If DFT_SCALE is set, the scaling is done after the transformation.
    // When DFT_COMPLEX_OUTPUT is set, the output is a complex matrix of the same size as input.
    cv::dft(img_fft, img_fft, cv::DFT_COMPLEX_OUTPUT); // cv::DFT_SCALE | 

    // # zero-out the center of the FFT shift (i.e., remove low
    // # frequencies), apply the inverse shift such that the DC
    // # component once again becomes the top-left, and then apply
    // # the inverse FFT

    // rearrange the quadrants of the result, so that the origin (zero, zero) corresponds with the image center.
    int cx = img_gray.cols / 2;
    int cy = img_gray.rows / 2;

    // center low frequencies in the middle
    // by shuffling the quadrants.
    cv::Mat q0(img_fft, cv::Rect(0, 0, cx, cy));    // Top-Left - Create a ROI per quadrant
    cv::Mat q1(img_fft, cv::Rect(cx, 0, cx, cy));   // Top-Right
    cv::Mat q2(img_fft, cv::Rect(0, cy, cx, cy));   // Bottom-Left
    cv::Mat q3(img_fft, cv::Rect(cx, cy, cx, cy));  // Bottom-Right

    cv::Mat tmp;  // swap quadrants (Top-Left with Bottom-Right)
    q0.copyTo(tmp);
    q3.copyTo(q0);
    tmp.copyTo(q3);

    q1.copyTo(tmp);  // swap quadrant (Top-Right with Bottom-Left)
    q2.copyTo(q1);
    tmp.copyTo(q2);

    if (debug_show) {
        std::vector<cv::Mat> planes;
        cv::Mat fft_mag;
        cv::split(img_fft, planes);                   // planes[0] = Re(DFT(I), planes[1] = Im(DFT(I))
        cv::magnitude(planes[0], planes[1], fft_mag);
        fft_mag += cv::Scalar::all(1); // switch to logarithmic scale
        cv::log(fft_mag, fft_mag);
        ShowColorMat("fft", fft_mag, true, 1, false, 0, 20);
    }

    // Block the low frequencies
    cv::mulSpectrums(img_fft, filter_mask, img_fft, 0); // multiply 2 spectrums

    if (debug_show) {
        std::vector<cv::Mat> planes;
        cv::Mat fft_mag;
        cv::split(img_fft, planes);                   // planes[0] = Re(DFT(I), planes[1] = Im(DFT(I))
        cv::magnitude(planes[0], planes[1], fft_mag);
        fft_mag += cv::Scalar::all(1); // switch to logarithmic scale
        cv::log(fft_mag, fft_mag);
        ShowColorMat("fft block low frequencies", fft_mag, true, 1, false, 0, 20);
    }

    // shuffle the quadrants to their original position
    cv::Mat p0(img_fft, cv::Rect(0, 0, cx, cy));    // Top-Left - Create a ROI per quadrant
    cv::Mat p1(img_fft, cv::Rect(cx, 0, cx, cy));   // Top-Right
    cv::Mat p2(img_fft, cv::Rect(0, cy, cx, cy));   // Bottom-Left
    cv::Mat p3(img_fft, cv::Rect(cx, cy, cx, cy));  // Bottom-Right

    p0.copyTo(tmp);
    p3.copyTo(p0);
    tmp.copyTo(p3);

    p1.copyTo(tmp);  // swap quadrant (Top-Right with Bottom-Left)
    p2.copyTo(p1);
    tmp.copyTo(p2);

    cv::dft(img_fft, img_fft, cv::DFT_SCALE | cv::DFT_INVERSE);

    std::vector<cv::Mat> complex_number;
    cv::Mat img_blur;
    cv::split(img_fft, complex_number); // planes[0] = Re(DFT(I), planes[1] = Im(DFT(I))
    magnitude(complex_number[0], complex_number[1], img_blur); // abs of complex number


    double min_val, max_val;
    cv::minMaxLoc(img_blur, &min_val, &max_val);

    if (max_val <= 0.f) {
        blur_value = 100;
        // black image
        return true;
    }

    cv::log(img_blur, img_blur);
    blur_value = cv::mean(img_blur)[0] * 20.f;

    return blur_value < thresh;
}

int main(int argc, char *argv[]) {
    float blur_thresh = 14; 
    float filter_mask_radius = 50;
    cv::Size mask_size = cv::Size(640, 360);

    cv::Mat img = cv::Mat::zeros(cv::Size(1280, 720), CV_8UC1);
    cv::Mat filter_mask = FilterMask(filter_mask_radius, mask_size);
    float blur_value;
    BlurDetection(img, filter_mask, blur_value, blur_thresh, true); 
    cv::imshow("input", img);
    std::cout<<blur_value<<std::endl;
    cv::waitKey(0);

    return 1;
}

#include <MNN/Interpreter.hpp>
#include <MNN/Matrix.h>
#include <MNN/ImageProcess.hpp>
#include <iostream>
#include <opencv2/opencv.hpp>
#include <opencv2/dnn/dnn.hpp>  //for cv::dnn::blobFromImage
#include <sys/time.h>

cv::Mat ShowMat(const cv::Mat& src) {
    double min;
    double max;
    cv::minMaxIdx(src, &min, &max);
    cv::Mat adjMap;

    float scale = 255 / (max - min);
    src.convertTo(adjMap, CV_8UC1, scale, -min * scale);

    cv::Mat falseColorsMap;
    cv::applyColorMap(adjMap, falseColorsMap, cv::COLORMAP_PINK);

    return falseColorsMap;
}

int main(int argc, char* argv[]) {
    std::string img_file = "/Users/luohanjie/Workspace/Vision/depth_estimation/MiDaS/input/squirrel_iphone_sample3.png";
    std::string model_file = "/Users/luohanjie/Workspace/Vision/my_slam/data/models/dpt_swin2_tiny_256/dpt_swin2_tiny_256.mnn";

    cv::Mat img = cv::imread(img_file);
    if (img.empty()) {
        std::cout << "Can not load image: " << img_file << std::endl;
        return 0;
    }

    int width_ori = img.cols;
    int height_ori = img.rows;

    // Interpreter是模型数据的持有者；Session通过Interpreter创建，是推理数据的持有者。多个推理可以共用同一个模型，即，多个Session可以共用一个Interpreter。
    // 在创建完Session，且不再创建Session或更新训练模型数据时，Interpreter可以通过releaseModel函数释放模型数据，以节省内存。
    std::shared_ptr<MNN::Interpreter> net(MNN::Interpreter::createFromFile(model_file.c_str()), MNN::Interpreter::destroy);
    if (net == NULL) {
        std::cout << "Can not load model: " << model_file << std::endl;
        return 0;
    }

    // 函数返回的Session实例是由Interpreter管理，随着Interpreter销毁而释放，一般不需要关注。也可以在不再需要时，调用Interpreter::releaseSession释放，减少内存占用。
    // 创建Session 一般而言需要较长耗时，而Session在多次推理过程中可以重复使用，建议只创建一次多次使用。
    MNN::ScheduleConfig session_config;
    session_config.type = MNN_FORWARD_AUTO;

    // memory、power、precision分别为内存、功耗和精度偏好。支持这些选项的后端会在执行时做出相应调整；若不支持，则忽略选项。
    // 示例： 后端 OpenCL precision 为 Low 时，使用 fp16 存储与计算，计算结果与CPU计算结果有少量误差，实时性最好；precision 为 Normal 时，使用 fp16存储，计算时将fp16转为fp32计算，计算结果与CPU计算结果相近，实时性也较好；precision 为 High 时，使用 fp32 存储与计算，实时性下降，但与CPU计算结果保持一致。
    // 后端 CPU precision 为 Low 时，根据设备情况开启 FP16 计算 precision 为 Low_BF16 时，根据设备情况开启 BF16 计算
    // BackendConfig bnconfig;
    // bnconfig.precision = BackendConfig::Precision_Low;
    // config.backendConfig = &bnconfig;
    MNN::Session* session = net->createSession(session_config);

    // 获取输入/出tensor
    MNN::Tensor* input = net->getSessionInput(session, "input.1");
    MNN::Tensor* output = net->getSessionOutput(session, "3335");

    // NCHW
    std::vector<int> input_dims = input->shape();
    int input_n = input_dims[0];
    int input_c = input_dims[1];
    int input_h = input_dims[2];
    int input_w = input_dims[3];
    std::cout << "Model input_n: "<<input_n<<", input_c: " << input_c<<", input_h: " << input_h << ", input_w: " << input_w << std::endl;

    // CHW
    std::vector<int> output_dims = output->shape();
    int output_c = output_dims[0];
    int output_h = output_dims[1];
    int output_w = output_dims[2];
    std::cout << "Model output_c: "<<output_c<<", output_h: " << output_h << ", output_w: " << output_w << std::endl;

    // x = (x / 255 - mean) / std
    // opencv: x = alpha * x + beta = (x / 255 - mean) / std = x / (255 * std) - mean / std
    // so:  alpha = 1 / (255 * std); beta = - mean / std
    float mean = 0.5f;
    float std = 0.5f;

    // N代表数量， C代表channel，H代表高度，W代表宽度。
    // NCHW其实代表的是[W H C N]， 第一个元素是000，第二个元素是沿着w方向的，即001，这样下去002 003，再接着呢就是沿着H方向，即004 005 006 007…这样到019后，沿C方向，轮到了020，之后021 022 …一直到319，然后再沿N方向。
    // NHWC代表的是[C W H N]， 第一个元素是000，第二个沿C方向，即020，040, 060…一直到300，之后沿W方向，001 021 041 061…301…到了303后，沿H方向，即004 024 …304.。最后到了319，变成N方向，320,340…
    // 当在不同的硬件加速的情况下，选用的类型不同，在intel GPU加速的情况下，因为GPU对于图像的处理比较多，希望在访问同一个channel的像素是连续的，一般存储选用NCHW，这样在做CNN的时候，在访问内存的时候就是连续的了，比较方便；
    // 所以在深度学习的时候，推理的前处理，一般都是将RGB或BGR图像进行转变为NCHW的格式；通常我们用opencv读取图像是NHWC的格式，需要进行通道分离，因为网路是一个通道一个通道的对图像做卷积，提取feature，所以NCHW更适合CNN。
    // https://blog.csdn.net/u010368556/article/details/105423260
    // caffe: NCHW;
    // pytorch: NCHW;
    // mxnet: NCHW;
    // 海思bgr: NCHW; 
    // NCNN: CHW
    // tensorflow: NHWC
    // opencv: NHWC
    // 瑞芯微rknn: NHWC
    // scipy.misc: NHW
    // https://www.cnblogs.com/yongy1030/p/11728103.html
    //convert NHWC to NCHW
    cv::Mat img_nchw;
    cv::dnn::blobFromImage(img, img_nchw,  1 / (255 * std), cv::Size(input_w, input_h), - mean / std, true);   //convert HWC to NCHW

    MNN::Tensor* tensor_nchw = new MNN::Tensor(input, MNN::Tensor::CAFFE);
    MNN::Tensor* tensor_depth = new MNN::Tensor(output, MNN::Tensor::CAFFE);

    memcpy(tensor_nchw->host<float>(), img_nchw.data, tensor_nchw->size());

    input->copyFromHostTensor(tensor_nchw);
    
    net->runSession(session);

    output->copyToHostTensor(tensor_depth);

    cv::Mat img_depth(output_h, output_w, CV_32FC1);   //difine opencv out img
    memcpy(img_depth.data, tensor_depth->host<float>(), tensor_depth->size());     //copy to output_img

    cv::resize(img_depth, img_depth, cv::Size(width_ori, height_ori));

    cv::Mat img_show = ShowMat(img_depth);
    cv::imshow("img_depth", img_show);
    cv::waitKey(0);

    delete tensor_nchw;
    delete tensor_depth;

    return 1;
}