1. warpaffine 介绍

放射变换(warpaffine)，主要解决图像的缩放和平移来处理目标检测中常见的预处理行为。

比如有一张猫的图片，对图片进行letterbox变换，将图片进行缩放，并且长边对齐，左右填充，这个时候用warpaffine就合适进行变换实现。

warpaffine说明

• 1. warpaffine是对图像做平移缩放旋转变换进行综合统一描述的方法
• 2. warpaffine也是一个很容易实现cuda并行加速的算法
• 3. 在深度学习领域通常需要做预处理，比如CopyMakeBorder(填充边)，RGB-BGR，减均值除以标准差，BGRBGRBGR->BBBGGGRRR
• 4. 如果使用cuda进行并行加速实现，那么可以对整个预处理都进行统一，并且性能非常好
     
     上图红框就是我们常见的图像预处理的代码，它需要做resize,然后复制填充边框，然后颜色空间BGR2RGB，然后再减均值除标准差，然后transpose变为tensor(batch,c,w,h) ,这5个步骤是非常耗时的，我们可以通过一个warpaffine变换，即可以cuda加速又可以并行一次性把5个步骤全部处理好，这样是我们利用warpaffine变换的主要目的
• 5. 由于warpaffine是标准的矩阵映射坐标，并且可逆，所以逆变换就是其变换矩阵的逆矩阵
• 6. 对于只有缩放和平移的变换矩阵，其有效自由度为3
     
     只有平移和算法，有效自由度就是scale,ox,oy,这3个值就可以了

2. 代码讲解

#include <cuda_runtime.h>
#include <opencv2/opencv.hpp>
#include <stdio.h>

using namespace cv;

#define min(a, b)  ((a) < (b) ? (a) : (b))
#define checkRuntime(op)  __check_cuda_runtime((op), #op, __FILE__, __LINE__)

bool __check_cuda_runtime(cudaError_t code, const char* op, const char* file, int line){
    if(code != cudaSuccess){
        const char* err_name = cudaGetErrorName(code);    
        const char* err_message = cudaGetErrorString(code);  
        printf("runtime error %s:%d  %s failed. \n  code = %s, message = %s\n", file, line, op, err_name, err_message);   
        return false;
    }
    return true;
}

void warp_affine_bilinear( // 声明
    uint8_t* src, int src_line_size, int src_width, int src_height, 
    uint8_t* dst, int dst_line_size, int dst_width, int dst_height, 
	uint8_t fill_value
);

Mat warpaffine_to_center_align(const Mat& image, const Size& size){  
    /*
       建议先阅读代码，若有疑问，可点击抖音短视频进行辅助讲解(建议1.5倍速观看)
            思路讲解：https://v.douyin.com/NhrNnVm/
            代码讲解: https://v.douyin.com/NhMv4nr/
    */        

    Mat output(size, CV_8UC3);
    uint8_t* psrc_device = nullptr;
    uint8_t* pdst_device = nullptr;
    size_t src_size = image.cols * image.rows * 3;
    size_t dst_size = size.width * size.height * 3;

    checkRuntime(cudaMalloc(&psrc_device, src_size)); // 在GPU上开辟两块空间
    checkRuntime(cudaMalloc(&pdst_device, dst_size));
    checkRuntime(cudaMemcpy(psrc_device, image.data, src_size, cudaMemcpyHostToDevice)); // 搬运数据到GPU上
    
    warp_affine_bilinear(
        psrc_device, image.cols * 3, image.cols, image.rows,
        pdst_device, size.width * 3, size.width, size.height,
        114
    );

    // 检查核函数执行是否存在错误
    checkRuntime(cudaPeekAtLastError());
    checkRuntime(cudaMemcpy(output.data, pdst_device, dst_size, cudaMemcpyDeviceToHost)); // 将预处理完的数据搬运回来
    checkRuntime(cudaFree(psrc_device));
    checkRuntime(cudaFree(pdst_device));
    return output;
}

int main(){ 
    /*
    若有疑问，可点击抖音短视频辅助讲解(建议1.5倍速观看) 
        https://v.douyin.com/NhMrb2A/
     */
    // int device_count = 1;
    // checkRuntime(cudaGetDeviceCount(&device_count));

    Mat image = imread("yq.jpg");
    Mat output = warpaffine_to_center_align(image, Size(640, 640));
    imwrite("output.jpg", output);
    printf("Done. save to output.jpg\n");
    return 0;
}

主要通过warpaffine_to_center_align函数实现warpaffine变换。

• 将数据搬到GPU上，开辟两块空间用来存放psrc_device,pdst_device数据

Mat output(size, CV_8UC3);
uint8_t* psrc_device = nullptr;
uint8_t* pdst_device = nullptr;
size_t src_size = image.cols * image.rows * 3;
size_t dst_size = size.width * size.height * 3;

checkRuntime(cudaMalloc(&psrc_device, src_size)); // 在GPU上开辟两块空间
checkRuntime(cudaMalloc(&pdst_device, dst_size));
checkRuntime(cudaMemcpy(psrc_device, image.data, src_size, cudaMemcpyHostToDevice)); // 搬运数据到GPU上

• 将psrc_device,pdst_device数据传入函数warp_affine_bilinear来处理

void warp_affine_bilinear(
    /* 
    建议先阅读代码，若有疑问，可点击抖音短视频进行辅助讲解(建议1.5倍速观看)
        - https://v.douyin.com/Nhre7fV/
     */
    uint8_t* src, int src_line_size, int src_width, int src_height, 
    uint8_t* dst, int dst_line_size, int dst_width, int dst_height, 
	uint8_t fill_value
){
    dim3 block_size(32, 32); // blocksize最大就是1024，这里用2d来看更好理解
    dim3 grid_size((dst_width + 31) / 32, (dst_height + 31) / 32);
    AffineMatrix affine;
    affine.compute(Size(src_width, src_height), Size(dst_width, dst_height));

    warp_affine_bilinear_kernel<<<grid_size, block_size, 0, nullptr>>>(
        src, src_line_size, src_width, src_height,
        dst, dst_line_size, dst_width, dst_height,
        fill_value, affine
    );
}

1. 函数里面调用了核函数warp_affine_bilinear_kernel, 定义了2d layout的grid_size和block_size, blocksize最大就是1024，32x32=1024刚好到了最大限值。grid_size定义为 dim3 grid_size((dst_width + 31) / 32, (dst_height + 31) / 32);
2. 定义了affine变换矩阵，并计算出变换矩阵，affine变换矩阵后续会转门写博客讲解
3. 调用warp_affine_bilinear_kernel核函数，把通过grid_size,block_size定义好的线程参数传进去。启动的线程数本质上就是dst图像的宽*高,只不过线程数由于向上取整会多一些，

 warp_affine_bilinear_kernel<<<grid_size, block_size, 0, nullptr>>>(
        src, src_line_size, src_width, src_height,
        dst, dst_line_size, dst_width, dst_height,
        fill_value, affine
    );

4. 核函数的代码

__global__ void warp_affine_bilinear_kernel(
    uint8_t* src, int src_line_size, int src_width, int src_height, 
    uint8_t* dst, int dst_line_size, int dst_width, int dst_height, 
	uint8_t fill_value, AffineMatrix matrix
){
    /* 
    建议先阅读代码，若有疑问，可点击抖音短视频进行辅助讲解(建议1.5倍速观看)
        - https://v.douyin.com/Nhr4vTF/
     */
    
    int dx = blockDim.x * blockIdx.x + threadIdx.x; 
    int dy = blockDim.y * blockIdx.y + threadIdx.y;
    if (dx >= dst_width || dy >= dst_height)  return;

    float c0 = fill_value, c1 = fill_value, c2 = fill_value;
    float src_x = 0; float src_y = 0;
    affine_project(matrix.d2i, dx, dy, &src_x, &src_y);

    
    /*
    建议先阅读代码，若有疑问，可点击抖音短视频进行辅助讲解(建议1.5倍速观看)
        - 双线性理论讲解：https://v.douyin.com/NhrH2tb/
        - 代码代码：https://v.douyin.com/NhrBqpc/ 
     */
    if(src_x < -1 || src_x >= src_width || src_y < -1 || src_y >= src_height){
        // out of range
        // src_x < -1时，其高位high_x < 0，超出范围
        // src_x >= -1时，其高位high_x >= 0，存在取值
    }else{
        int y_low = floorf(src_y);
        int x_low = floorf(src_x);
        int y_high = y_low + 1;
        int x_high = x_low + 1;

        uint8_t const_values[] = {fill_value, fill_value, fill_value};
        float ly    = src_y - y_low;
        float lx    = src_x - x_low;
        float hy    = 1 - ly;
        float hx    = 1 - lx;
        float w1    = hy * hx, w2 = hy * lx, w3 = ly * hx, w4 = ly * lx;
        uint8_t* v1 = const_values;
        uint8_t* v2 = const_values;
        uint8_t* v3 = const_values;
        uint8_t* v4 = const_values;
        if(y_low >= 0){
            if (x_low >= 0)
                v1 = src + y_low * src_line_size + x_low * 3;

            if (x_high < src_width)
                v2 = src + y_low * src_line_size + x_high * 3;
        }
        
        if(y_high < src_height){
            if (x_low >= 0)
                v3 = src + y_high * src_line_size + x_low * 3;

            if (x_high < src_width)
                v4 = src + y_high * src_line_size + x_high * 3;
        }
        
        c0 = floorf(w1 * v1[0] + w2 * v2[0] + w3 * v3[0] + w4 * v4[0] + 0.5f);
        c1 = floorf(w1 * v1[1] + w2 * v2[1] + w3 * v3[1] + w4 * v4[1] + 0.5f);
        c2 = floorf(w1 * v1[2] + w2 * v2[2] + w3 * v3[2] + w4 * v4[2] + 0.5f);
    }

    uint8_t* pdst = dst + dy * dst_line_size + dx * 3;
    pdst[0] = c0; pdst[1] = c1; pdst[2] = c2;
}

2D layout的线程，相当于2个1D，从它的定义可以看出，x方向和y方向的线程索引定义方式是一样的，只不过将x变为y,所以说是两个1D

int dx = blockDim.x * blockIdx.x + threadIdx.x; 
int dy = blockDim.y * blockIdx.y + threadIdx.y;

当x,y方向线程索引超过图片宽高，则直接就return返回

if (dx >= dst_width || dy >= dst_height)  return;

拿到x,y方向的线程索引dx,dy，其实就相当于目标图dst 的x,y坐标，然后通过矩阵映射affine_project,映射为输入图的x,y坐标，拿到输入的x,y坐标就可以做双线性插值。

float c0 = fill_value, c1 = fill_value, c2 = fill_value;
    float src_x = 0; float src_y = 0;
    affine_project(matrix.d2i, dx, dy, &src_x, &src_y);

    
    /*
    建议先阅读代码，若有疑问，可点击抖音短视频进行辅助讲解(建议1.5倍速观看)
        - 双线性理论讲解：https://v.douyin.com/NhrH2tb/
        - 代码代码：https://v.douyin.com/NhrBqpc/ 
     */
    if(src_x < -1 || src_x >= src_width || src_y < -1 || src_y >= src_height){
        // out of range
        // src_x < -1时，其高位high_x < 0，超出范围
        // src_x >= -1时，其高位high_x >= 0，存在取值
    }else{
        int y_low = floorf(src_y);
        int x_low = floorf(src_x);
        int y_high = y_low + 1;
        int x_high = x_low + 1;

        uint8_t const_values[] = {fill_value, fill_value, fill_value};
        float ly    = src_y - y_low;
        float lx    = src_x - x_low;
        float hy    = 1 - ly;
        float hx    = 1 - lx;
        float w1    = hy * hx, w2 = hy * lx, w3 = ly * hx, w4 = ly * lx;
        uint8_t* v1 = const_values;
        uint8_t* v2 = const_values;
        uint8_t* v3 = const_values;
        uint8_t* v4 = const_values;
        if(y_low >= 0){
            if (x_low >= 0)
                v1 = src + y_low * src_line_size + x_low * 3;

            if (x_high < src_width)
                v2 = src + y_low * src_line_size + x_high * 3;
        }
        
        if(y_high < src_height){
            if (x_low >= 0)
                v3 = src + y_high * src_line_size + x_low * 3;

            if (x_high < src_width)
                v4 = src + y_high * src_line_size + x_high * 3;
        }

算完后，就得到了3个像素值c0,c1,c2,就是对应dx,dy处应该给定的像素值,然后再填充到目标图pdst中

        c0 = floorf(w1 * v1[0] + w2 * v2[0] + w3 * v3[0] + w4 * v4[0] + 0.5f);
        c1 = floorf(w1 * v1[1] + w2 * v2[1] + w3 * v3[1] + w4 * v4[1] + 0.5f);
        c2 = floorf(w1 * v1[2] + w2 * v2[2] + w3 * v3[2] + w4 * v4[2] + 0.5f);

由于图像有多大，我们就启动了多大的线程，所以线程时针对每个像素处理的。如果这里想希望实现BGR2RGB，你会发现你只需要做顺序交换就行了。同样减均值除标准差也可以基于每个线程对每个像素进行操作：pixel -mean/std

uint8_t* pdst = dst + dy * dst_line_size + dx * 3;
pdst[0] = c0; pdst[1] = c1; pdst[2] = c2;

// BGR -> RGB
pdst[2]=c0;pdst[1]=c1;pdst[0]=c2;

代码下载：

链接：https://pan.baidu.com/s/1wqPF40dHd7yzcXjLFYY5yw?pwd=bi5n
提取码：bi5n