上个月做完一个图像处理的IP核设计，由于涉及到上采样的超分辨率算法，就看了一下AMD 开源的超分算法——FSR（FidelityFX Super Resolution），并打算用C源代码实现（without lib）FSR中的上采样（EASU）部分。
FSR包括两个部分（两个Pass）：上采样（EASU）+ 锐化（RCAS），而我打算主要是实现EASU部分，至于FSR的数学逻辑和算法介绍，可以参考以下文章，非常详细！！！！
详细剖析 AMD FSR 算法 - 知乎https://zhuanlan.zhihu.com/p/401030221
温馨提示：在看我的代码之前，希望小可爱们可以提前学习下图像处理中最基本的算法，比如灰度转换，边缘检测，高斯滤波等等，还有就是C语言中的基本知识哦。废话不多说，直接上代码 ，go ~~~~~~
2.代码实现 C Model:#include支持函数：bmp.c#include #include #include #include "bmp.h" //全局变量 float FXe,FYe,FXd,FYd,FXh,FYh,FXi,FYi,Fe,Fd,Fh,Fi,We,Wd,Wh,Wi,W; float dx,dy,feature,Sy,Sx,lobs,clps,n_x,n_y; float clamp(float val,float min, float max) { if(val<=min) val = min; else if(val>=max) val = max; else val = val; return val; } // (25/16 * (2/5 * x^2 - 1)^2 - (25/16 - 1)) * (w * x^2 - 1)^2 // |_______________________________________| |_______________| // base window float Lanczos2(float off_x,float off_y, float Dir_x, float Dir_y, float Len_x, float Len_y, float Lob, float Clp) { float L_w,dxx,dyy,x_number; dxx = (off_x * Dir_x) + (off_y * Dir_y); dyy = (off_x * -1 * Dir_y) + (off_y * Dir_x); dxx = dxx * Len_x; dyy = dyy * Len_y; x_number = dxx * dxx + dyy * dyy; //做distance^2 x_number = fmin(x_number,Clp); //裁剪 float WB = 0.4f * x_number - 1.0f; float WA = Lob * x_number - 1.0f; WB *= WB; WA *= WA; L_w = 1.5625f * WB + (-0.5625f); L_w = L_w * WA; return L_w; } void distance(uint32_t x, uint32_t y, float i, float j) { n_x = x - i; //输出off_x值 n_y = -(j - y); } float gradient(BMPImage *src, uint32_t x, uint32_t y, float i, float j) { float Weight_e = (x+1-i)*(y+1-j); float Weight_d = (i-x)*(y+1-j); float Weight_h = (i-x)*(j-y); float Weight_i = (x+1-i)*(j-y); float DirR,Stretch; uint8_t Zero; Pixel *in1_a = bmp_pixel_at(src , x , y-1); Pixel *in1_b = bmp_pixel_at(src , x+1 , y-1); Pixel *in1_f = bmp_pixel_at(src , x-1 , y); Pixel *in1_e = bmp_pixel_at(src , x , y); Pixel *in1_d = bmp_pixel_at(src , x+1 , y); Pixel *in1_c = bmp_pixel_at(src , x+2 , y); Pixel *in1_j = bmp_pixel_at(src , x-1 , y+1); Pixel *in1_i = bmp_pixel_at(src , x , y+1); Pixel *in1_h = bmp_pixel_at(src , x+1 , y+1); Pixel *in1_g = bmp_pixel_at(src , x+2 , y+1); Pixel *in1_k = bmp_pixel_at(src , x , y+2); Pixel *in1_n = bmp_pixel_at(src , x+1 , y+2); FXe = abs(in1_f->r - in1_d->r) / fmax(abs(in1_f->r - in1_e->r),abs(in1_d->r - in1_e->r)); FYe = abs(in1_a->r - in1_i->r) / fmax(abs(in1_a->r - in1_e->r),abs(in1_i->r - in1_e->r)); dx += (in1_d->r - in1_f->r) * Weight_e; dy += (in1_i->r - in1_a->r) * Weight_e; FXe = clamp(FXe,0.0f,1.0f); FYe = clamp(FYe,0.0f,1.0f); FXe *= FXe; FYe *= FYe; feature += (FXe + FYe) * Weight_e; FXd = abs(in1_e->r - in1_c->r) / fmax(abs(in1_e->r - in1_d->r),abs(in1_c->r - in1_d->r)); FYd = abs(in1_b->r - in1_h->r) / fmax(abs(in1_b->r - in1_d->r),abs(in1_h->r - in1_d->r)); dx += (in1_c->r - in1_e->r) * Weight_d; dy += (in1_h->r - in1_b->r) * Weight_d; FXd = clamp(FXd,0.0f,1.0f); FYd = clamp(FYd,0.0f,1.0f); FXd *= FXd; FYd *= FYd; feature += (FXd + FYd) * Weight_d; FXh = abs(in1_i->r - in1_g->r) / fmax(abs(in1_i->r - in1_h->r),abs(in1_g->r - in1_h->r)); FYh = abs(in1_d->r - in1_n->r) / fmax(abs(in1_d->r - in1_h->r),abs(in1_n->r - in1_h->r)); dx += (in1_g->r - in1_i->r) * Weight_h; dy += (in1_n->r - in1_d->r) * Weight_h; FXh = clamp(FXh,0.0f,1.0f); FYh = clamp(FYh,0.0f,1.0f); FXh *= FXh; FYh *= FYh; feature += (FXh + FXh) * Weight_h; FXi = abs(in1_j->r - in1_h->r) / fmax(abs(in1_j->r - in1_i->r),abs(in1_h->r - in1_i->r)); FYi = abs(in1_e->r - in1_k->r) / fmax(abs(in1_e->r - in1_i->r),abs(in1_k->r - in1_i->r)); dx += (in1_h->r - in1_j->r) * Weight_i; dy += (in1_k->r - in1_e->r) * Weight_i; FXi = clamp(FXi,0.0f,1.0f); FYi = clamp(FYi,0.0f,1.0f); FXi *= FXi; FYi *= FYi; feature += (FXi + FXi) * Weight_i; //梯度归一化 DirR = dx * dx + dy * dy; Zero = DirR < (1.0f/32768.0) ? 1 : 0; DirR = 1.0f / sqrt(DirR); DirR = Zero ? 1.0f : DirR; dx = Zero ? 1.0f : dx; dy = Zero ? 1.0f : dy; dx *= DirR; //cos(theta),输出值 dy *= DirR; //sin(theta)输出值 //加入空间调整 feature = feature * 0.5f; feature *= feature; //输出feature值 Stretch = 1.0f / fmax(abs(dx),abs(dy)); Sx = 1.0 + (Stretch - 1.0) * feature; //输出调整空间 Sy= 1.0 - 0.5 * feature; lobs = 0.5 - 0.29 * feature; //0.5 - 0.25 * feature = lobs clps = 1.0f / lobs; //输出Lob和Clp值 } int QIT(BMPImage *img , uint32_t x ,uint32_t y) { double W_in_0 = (img->header.width_px +3 ) / 4 *4; double H_in_0 = img->header.height_px; uint8_t i , Q; Pixel *out_pixel = bmp_pixel_at(img , x , y); i = out_pixel->r; if(i==255) Q = 1; else Q = 0; return Q; } int boundray(uint32_t x, uint32_t y,uint32_t i,uint32_t j) { uint8_t CT; if((x>=7 && x<=i-8) && (y>=7 && y<=j-8)) CT = 1; else CT = 0; return CT; } BMPImage *FSR(BMPImage *img0, BMPImage *img1, BMPImage *img3, BMPImage *img2) //FSR { double W_in_0 = (img1->header.width_px +3 ) / 4 *4; double H_in_0 = img1->header.height_px; uint32_t i,j; BMPImage *dst = bmp_create(W_in_0, H_in_0); //创建目标图像 for(j = 0; j< H_in_0; j++) for(i = 0; i< W_in_0; i++) { uint8_t WT = boundray(i,j,W_in_0,H_in_0); if(WT==1) { uint8_t QT = QIT(img3 , i , j); //表示为平坦点 if(QT == 0) { Pixel *out_dst = bmp_pixel_at(dst , i , j); Pixel *out_img1 = bmp_pixel_at(img1 , i , j); out_dst-> r = out_img1 ->r; out_dst-> g = out_img1 ->g; out_dst-> b = out_img1 ->b; } //表示边缘点 else { Pixel *out_dst1 = bmp_pixel_at(dst , i , j); double F_W = i * 0.25f + 0.5f * (0.25f -1.0f); //浮点数中心像素(y,x) double F_H = j * 0.25f + 0.5f * (0.25f -1.0f); uint32_t x = floor(F_W); //浮点像素向下取整 uint32_t y = floor(F_H); float La,Lb,Lf,Le,Ld,Lc,Lj,Li,Lh,Lg,Lk,Ln,L_sum; dx = 0;dy = 0;feature = 0; gradient(img2, x, y, F_W, F_H); //一个 distance(x, y-1, F_W, F_H); //共12组 La = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x+1, y-1, F_W, F_H); Lb = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x-1, y, F_W, F_H); Lf = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x, y, F_W, F_H); Le = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x+1, y, F_W, F_H); Ld = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x+2, y, F_W, F_H); Lc = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x-1, y+1, F_W, F_H); Lj = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x, y+1, F_W, F_H); Li = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x+1, y+1, F_W, F_H); Lh = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x+2, y+1, F_W, F_H); Lg = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x, y+2, F_W, F_H); Lk = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); distance(x+1, y+2, F_W, F_H); Ln = Lanczos2(n_x, n_y, dx, dy, Sx, Sy, lobs, clps); L_sum = La+Lb+Lf+Le+Ld+Lc+Lj+Li+Lh+Lg+Lk+Ln; Pixel *in2_a = bmp_pixel_at(img0 , x , y-1); Pixel *in2_b = bmp_pixel_at(img0 , x+1 , y-1); Pixel *in2_f = bmp_pixel_at(img0 , x-1 , y); Pixel *in2_e = bmp_pixel_at(img0 , x , y); Pixel *in2_d = bmp_pixel_at(img0 , x+1 , y); Pixel *in2_c = bmp_pixel_at(img0 , x+2 , y); Pixel *in2_j = bmp_pixel_at(img0 , x-1 , y+1); Pixel *in2_i = bmp_pixel_at(img0 , x , y+1); Pixel *in2_h = bmp_pixel_at(img0 , x+1 , y+1); Pixel *in2_g = bmp_pixel_at(img0 , x+2 , y+1); Pixel *in2_k = bmp_pixel_at(img0 , x , y+2); Pixel *in2_n = bmp_pixel_at(img0 , x+1 , y+2); uint8_t A_r = fmax(in2_e->r , fmax(in2_d->r , fmax(in2_h->r , in2_i->r))); uint8_t A_g = fmax(in2_e->g , fmax(in2_d->g , fmax(in2_h->g , in2_i->g))); uint8_t A_b = fmax(in2_e->g , fmax(in2_d->g , fmax(in2_h->g , in2_i->g))); uint8_t B_r = fmin(in2_e->r , fmin(in2_d->r , fmin(in2_h->r , in2_i->r))); uint8_t B_g = fmin(in2_e->g , fmin(in2_d->g , fmin(in2_h->g , in2_i->g))); uint8_t B_b = fmin(in2_e->g , fmin(in2_d->g , fmin(in2_h->g , in2_i->g))); uint8_t rr = (in2_a->r * La + in2_b->r * Lb + in2_f->r * Lf + in2_e->r * Le + in2_d->r * Ld + in2_c->r * Lc + in2_j->r * Lj + in2_i->r * Li + in2_h->r * Lh + in2_g->r * Lg + in2_k->r * Lk + in2_n->r * Ln) / L_sum; uint8_t gg = (in2_a->g * La + in2_b->g * Lb + in2_f->g * Lf + in2_e->g * Le + in2_d->g * Ld + in2_c->g * Lc + in2_j->g * Lj + in2_i->g * Li + in2_h->g * Lh + in2_g->g * Lg + in2_k->g * Lk + in2_n->g * Ln) / L_sum; uint8_t bb = (in2_a->b * La + in2_b->b * Lb + in2_f->b * Lf + in2_e->b * Le + in2_d->b * Ld + in2_c->b * Lc + in2_j->b * Lj + in2_i->b * Li + in2_h->b * Lh + in2_g->b * Lg + in2_k->b * Lk + in2_n->b * Ln) / L_sum; out_dst1->r = fmin(A_r, fmax(B_r , rr)); out_dst1->g = fmin(A_g, fmax(B_g , gg)); out_dst1->b = fmin(A_b, fmax(B_b , bb)); } } else { Pixel *out1_dst = bmp_pixel_at(dst , i , j); Pixel *out1_img1 = bmp_pixel_at(img1 , i , j); out1_dst->r = out1_img1->r; out1_dst->g = out1_img1->g; out1_dst->b = out1_img1->b; } } return dst; }
#include#include #include #include "bmp.h" BMPImage *bmp_read(const char *filename) { FILE *fp = fopen(filename, "r"); if (!fp) { printf("file %s open failed!n", filename); goto OPEN_FAIL; } BMPImage *img = malloc(sizeof(BMPImage)); if (!img) { printf("image malloc fail!n"); goto IMG_MALLOC_FAIL; } fread(&img->header, 1, sizeof(BMPHeader), fp); if (img->header.type != 0x4D42) { printf("illegal bmp file %s", filename); goto NOT_BMP; } if (!img->header.image_size_bytes) { img->header.image_size_bytes = img->header.size - img->header.offset; } Pixel *pixel_data = malloc(img->header.image_size_bytes); if (!pixel_data) { printf("pixel data malloc fail!n"); goto PIX_MALLOC_FAIL; } rewind(fp); fread(pixel_data, 1, img->header.offset, fp); fread(pixel_data, 1, img->header.image_size_bytes, fp); img->data = pixel_data; fclose(fp); return img; PIX_MALLOC_FAIL: NOT_BMP: free(img); IMG_MALLOC_FAIL: fclose(fp); OPEN_FAIL: return 0; } BMPImage *bmp_create(uint32_t width, uint32_t height) { BMPImage *img = malloc(sizeof(BMPImage)); if (!img) { printf("img create malloc fail!n"); return 0; } memset(img, 0, sizeof(BMPImage)); img->data = malloc(width * height * 3); if (!img->data) { printf("img create data malloc fail!n"); return 0; } BMPHeader *h = &img->header; h->type = 0x4D42; h->size = width * height * 3 + 54; h->offset = 54; h->dib_header_size = 40; h->width_px = width; h->height_px = height; h->num_planes = 1; h->bits_per_pixel = 24; h->compression = 0; h->image_size_bytes = width * height * 3; return img; } int32_t bmp_write(BMPImage *img, const char *filename) { FILE *fp = fopen(filename, "wb+"); if (!fp) { printf("open/create file %s failed!n", filename); return -1; } fwrite(&img->header, 1, sizeof(BMPHeader), fp); fwrite(img->data, 1, img->header.image_size_bytes, fp); return 0; } Pixel *bmp_pixel_at(BMPImage *img, uint32_t x, uint32_t y) { if (x >= img->header.width_px || y >= img->header.height_px) { return 0; } return (Pixel *)((uint8_t *)img->data + (x + y * img->header.width_px) * img->header.bits_per_pixel / 8); } void bmp_print_header(BMPImage *img) { BMPHeader *h = &img->header; printf("type=%xn", h->type); printf("size=%un", h->size); printf("reserved1=%un", h->reserved1); printf("reserved2=%un", h->reserved2); printf("offset=%un", h->offset); printf("dib_header_size=%un", h->dib_header_size); printf("width_px=%dn", h->width_px); printf("height_px=%dn", h->height_px); printf("num_planes=%un", h->num_planes); printf("bits_per_pixel=%un", h->bits_per_pixel); printf("compression=%un", h->compression); printf("image_size_bytes=%un", h->image_size_bytes); printf("x_resolution_ppm=%dn", h->x_resolution_ppm); printf("y_resolution_ppm=%dn", h->y_resolution_ppm); printf("num_colors=%un", h->num_colors); printf("important_colors=%un", h->important_colors); }
支持函数：gray.c
#include#include #include "bmp.h" BMPImage *bmp_gray(BMPImage *img) { uint32_t w_i = (img->header.width_px + 3) / 4 * 4; uint32_t h_i = img->header.height_px; BMPImage *dst = bmp_create(w_i, h_i); uint32_t gray; // color channels uint32_t x, y; for (y = 0; y < h_i; y++) { for (x = 0; x < w_i; x++) { Pixel *p = bmp_pixel_at(img, x, y); gray = (p->r * 114 + p->g * 587 + p->b * 299) / 1000; Pixel *t = bmp_pixel_at(dst, x, y); t->r = (uint8_t)gray; t->g = (uint8_t)gray; t->b = (uint8_t)gray; } } return dst; }
支持函数：sobel.c
#include#include #include #include #include "bmp.h" #include "guass.h" #define SOBEL_THRESHOLD 128 //阈值，待优化 #define ROW 2160 #define COL 3840 uint32_t image_data[COL + 2][ROW + 2]; uint32_t Gxy[COL][ROW]; uint32_t Gx[COL][ROW]; uint32_t Gy[COL][ROW]; uint32_t THRESHOLD[COL][ROW]; void store(BMPImage *img, BMPImage *rst) { uint8_t MAX_ = 254; uint8_t MIN_ = 1; int i; uint32_t w_i = (img->header.width_px + 3) / 4 * 4; uint32_t h_i = img->header.height_px; uint32_t x, y; for (y = 1; y < h_i + 1; y++) { for (x = 1; x < w_i + 1; x++) { Pixel *p = bmp_pixel_at(rst, x - 1, y - 1); if (p->r >= MAX_) image_data[x][y] = 255; else if (p->r <= MIN_) image_data[x][y] = 0; else image_data[x][y] = p->r; } } } //填充图像扩充边沿两边为零 void padding(BMPImage *img) { int i; for (i = 0; i < img->header.width_px; i++) { image_data[i][0] = 0; image_data[i][img->header.height_px + 2] = 0; } for (i = 0; i < img->header.height_px; i++) { image_data[0][i] = 0; image_data[img->header.width_px + 2][i] = 0; } } //卷积3*3核 //中心像素(row,col) int convolution(int kernel[3][3], int row, int col) { int i, j, sum = 0; for (i = -1; i < 2; i++) { for (j = -1; j < 2; j++) { sum += image_data[j + col][i + row] * kernel[i+1][j+1]; } } return sum; } //sobel 算子边沿检测 void sobel_edge_detector(BMPImage *img) { uint32_t i, j, gx, gy; int mx[3][3] = { {-1, 0, 1}, {-2, 0, 2}, {-1, 0, 1} }; int my[3][3] = { {-1, -2, -1}, {0, 0, 0}, {1, 2, 1} }; uint32_t w_i = (img->header.width_px + 3) / 4 * 4; uint32_t h_i = img->header.height_px; for (i = 1; i < img->header.height_px + 1; i++) { for (j = 1; j < img->header.width_px + 1; j++) { gx = convolution(mx, i, j); gy = convolution(my, i, j); Gxy[j - 1][i - 1] = (uint8_t)sqrt(gx*gx + gy * gy); Gx[j - 1][i - 1] = (uint8_t)gx; Gy[j - 1][i - 1] = (uint8_t)gy; } } } //标准化 void min_max_normalization(BMPImage *img) { int min = 1000000, max = 0, i, j; uint32_t w_i = (img->header.width_px + 3) / 4 * 4; uint32_t h_i = img->header.height_px; for (i = 0; i < h_i; i++) { for (j = 0; j < w_i; j++) { if (Gxy[j][i] < min) { min = Gxy[j][i]; } else if (Gxy[j][i] > max) { max = Gxy[j][i]; } else { min = min; max = max; } } } for (i = 0; i < h_i; i++) { for (j = 0; j < w_i; j++) { double ratio = (double)(Gxy[j][i] - min) / (max - min); Gxy[j][i] = ratio * 255; } } } //二值化 void Image_Binarization(BMPImage *img) { int i, j; uint32_t w_i = (img->header.width_px + 3) / 4 * 4; uint32_t h_i = img->header.height_px; for (i = 0; i < h_i; i++) { for (j = 0; j < w_i; j++) { if (Gxy[j][i] > SOBEL_THRESHOLD) THRESHOLD[j][i] = 255; else THRESHOLD[j][i] = 0; } } } BMPImage *bmp_sobel(BMPImage *img,BMPImage *dst) { uint32_t w = (img->header.width_px + 3) / 4 * 4; uint32_t h = img->header.height_px; // BMPImage *dst0 = bmp_create(w, h); //创建dst图片 // BMPImage *dst = bmp_gray(img); //转灰度 // BMPImage *rst = GaussianFilter(img,dst); store(img, dst); //存储矩阵 padding(img); //边沿填充 sobel_edge_detector(img); //边缘检测 min_max_normalization(img); //标准化 Image_Binarization(img); //二值化 BMPImage *dst1 = bmp_create(w, h); //创建彩色边缘检测输出图片 int i, j; for (i = 0; i < h; i++) { for (j = 0; j < w; j++) { Pixel *s = bmp_pixel_at(dst1, j, i); s->r = (uint8_t)THRESHOLD[j][i]; s->g = (uint8_t)THRESHOLD[j][i]; s->b = (uint8_t)THRESHOLD[j][i]; } } return dst1; }
头文件应该问题不大，我就不放上去了，Main.c可以根据自己的想法搭建，我是在GCC上运行的，注意在gcc上运行这些文件，需要先产生 .o 文件，才能编译。
3.结果 由于是BMP位图格式，24位的，数据有点大，4K上传不了，所以我转成jpg格式来对比 低分辨率960x540 高分辨率3840x2160
