代码中集成了几种基于opencv的插值方法，自己写了一个自适应插值法，对比了效果。背景和原理介绍参考我另一篇博客：原理解读
一、MATLAB代码
1、2已经可以了。3、4还没完全调通,因为没太理解matlab里这个样条插值，3、4做出来的结果比原图还差。。。
这个是自适应的插值，与传统的插值方法不同，按理说应该会比直接opencv插效果要好，后面有对比。
代码处理是基于一个灰度的处理，如果是RGB图像的话可以调用三次再整合。
clear all;
X=imread('/home/yhy/图片/butterfly_GT.bmp');
f1=X;
f1=double(f1(:,:,1))/255;
[M,N]=size(f1);
h=[0.5,0.5];
g=filter(h,1,f1');
%g=g(1:2:N,:);
%g=filter(h,1,g');
%g=g(1:2:M,:);
%#SNR=input('SNR,DB');
%SNR=20
%gg=im2col(g,[(M)/2,N/2],'distinct');
%n_var=var(gg)/10^(SNR/10); #fangcha
%g=imnoise(g,'gaussian',0,n_var(:,1));
subplot(3,1,1);
imshow(g);
title(size(g));
#以上为图像缩小处理
#第一次竖向插值
%key0=input('press 1 for Bilinear,2 for Bicudic and 3 for Cubic Spline and 4 for cubic Moms:');
key0=2
#1：双线性；2：双立方；3：三次样条；4：三次样条OMOS
f=g;
[Mg,Ng]=size(g);
#if(key0==3)||(key0==4)
#fl=spline(0:(Ng-1),g); #spline(x,y)返回的结构体包含了三次样条插值所得的分段多项式式的相关信息
#fl=fl.coefs;
#f(:,:,1)=[reshape(fl(:,1),Mg,Ng-1),g(:,Ng)];
#fs=g;
#end;
a=-1/2;
[M,N]=size(g);
ff=zeros(M,N);
x=f(:,N-1:N);
x=rot90(x,2);
y=f(:,1);
f=[y,f,x];
L=256;
k=1.5;
s=1/2;
s=linspace(0,1,11);
l1=l2=l3=l4=1;
temp1=zeros(11,2*N);
h=[0.5,0.5];
for i=1:M
for c=1:11
for j=2:N+1
A=(abs(f(i,j+1)-f(i,j-1))-abs(f(i,j+2)-f(i,j)))/(L-1);
switch key0
case 1
temp(c,j-1)=l2*f(i,j)*(1-s(c))+l3*f(i,j+1)*s(c);
case 2
temp(c,j-1)=l1*f(i,j-1)*(a*s(c)^3-2*a*s(c)^2+a*s(c))+l2*f(i,j)*((a+2)*s(c)^3-(3+a)*s(c)^2+1)+l3*f(i,j+1)*(-(a+2)*s(c)^3+(2*a+3)*s(c)^2-a*s(c))+l4*f(i,j+2)*(-a*s(c)^3+a*s(c)^2);
case 3
temp(c,j-1)=l1*f(i,j-1)*((3+s(c))^3-4*(2+s(c))^3+6*(1+s(c))^3-4*s(c)^3)/6+l2*f(i,j)*((2+s(c))^3-4*(1+s(c))^3+6*s(c)^3)/6+l3*f(i,j+1)*((1+s(c))^3-4*s(c)^3)/6+l4*f(i,j+2)*s(c)^3/6 ;
case 4
temp(c,j-1)=l1*f(i,j-1)*((-1/6)*(1+s(c))^3+(1+s(c))^2+(-85/42)*(1+s(c))+(29/21))+l2*f(i,j)*(0.5*s(c)^3-s(c)^2+(1/14)*s(c)+13/21)+l3*f(i,j+1)*(0.5*(1-s(c))^3-(1-s(c))^2+(1/14)*s(c)+13/21)+l4*f(i,j+2)*((-1/6)*(2-s(c))^3+(2-s(c))^2-(85/42)*(2-s(c))+29/21);
end
end
end
for c=1:11
temp1(c,1:2:2*N)=f(i,2:N+1);
#if(key0==3)||(key0==4)
#temp1(c,1:2:2*N)=fs(i,1:N);
#end
end
temp1(1:11,2:2:2*N)=temp(1:11,1:N);
temp1=[temp1(:,1),temp1];
temp2=filter(h,1,temp1');
temp2=temp2';
temp3=temp2(:,1:2:2*N);
for c=1:11
error(c,:)=abs(f(i,2:N+1)-temp3(c,:));
#if(key0==3)||(key0==4)
#error(c,:)=abs(fs(i,1:N)-temp3(c,:));
#end
end
[min_val,min_idx]=min(error);
for j=1:N
ff(i,j)=temp(min_idx(j),j);
end
end
ff=ff(:,1:N);
fff(1:M,1:2:2*N)=f(1:M,2:N+1);
#if(key0==3)||(key0==4)
#fff(1:M,1:2:2*N)=fs(1:M,1:N);
#end
fff(1:M,2:2:2*N)=ff(1:M,1:N);
fff=(fff>=0).*fff;
f=fff';
subplot(3,1,2);
imshow(f);
title(size(f));
save data f;
#第二次横向
clear all
load data f;
%key0=input('press 1 for Bilinear,2 for Bicudic and 3 for Cubic Spline and 4 for cubic Moms:');
key0=2
#1：双线性；2：双立方；3：三次样条；4：三次样条OMOS
[M,N]=size(f);
#if(key0==3)||(key0==4)
#fs=f;
#fl=spline(0:M-1,f); #spline(x,y)返回的结构体包含了三次样条插值所得的分段多项式式的相关信息
#fl=fl.coefs;
#f(:,:)=[reshape(fl(:,1),M-1,N);fs(M,:)];
#end;
a=-1/2;
ff=zeros(M,N);
x=f(:,N-1:N);
x=rot90(x,2); #转向到横向
y=f(:,1);
f=[y,f,x];
L=256;
k=1.5;
k1=1.5;
s=1/2;
s=linspace(0,1,11);
l1=l2=l3=l4=1;
temp1=zeros(11,2*N);
h=[0.5,0.5];
for i=1:M
for c=1:11
for j=2:N+1
A=(abs(f(i,j+1)-f(i,j-1))-abs(f(i,j+2)-f(i,j)))/(L-1);
switch key0
case 1
temp(c,j-1)=l2*f(i,j)*(1-s(c))+l3*f(i,j+1)*s(c);
case 2
temp(c,j-1)=l1*f(i,j-1)*(a*s(c)^3-2*a*s(c)^2+a*s(c))+l2*f(i,j)*((a+2)*s(c)^3-(3+a)*s(c)^2+1)+l3*f(i,j+1)*(-(a+2)*s(c)^3+(2*a+3)*s(c)^2-a*s(c))+l4*f(i,j+2)*(-a*s(c)^3+a*s(c)^2);
case 3
temp(c,j-1)=l1*f(i,j-1)*((3+s(c))^3-4*(2+s(c))^3+6*(1+s(c))^3-4*s(c)^3)/6+l2*f(i,j)*((2+s(c))^3-4*(1+s(c))^3+6*s(c)^3)/6+l3*f(i,j+1)*((1+s(c))^3-4*s(c)^3)/6+l4*f(i,j+2)*s(c)^3/6 ;
case 4
temp(c,j-1)=l1*f(i,j-1)*((-1/6)*(1+s(c))^3+(1+s(c))^2+(-85/42)*(1+s(c))+(29/21))+l2*f(i,j)*(0.5*s(c)^3-s(c)^2+(1/14)*s(c)+13/21)+l3*f(i,j+1)*(0.5*(1-s(c))^3-(1-s(c))^2+(1/14)*s(c)+13/21)+l4*f(i,j+2)*((-1/6)*(2-s(c))^3+(2-s(c))^2-(85/42)*(2-s(c))+29/21);
end
end
end
for c=1:11
temp1(c,1:2:2*N)=f(i,2:N+1);
#if(key0==3)||(key0==4)
#temp1(c,1:2:2*N)=fs(i,1:N);
#end
end
temp1(1:11,2:2:2*N)=temp(1:11,1:N);
temp1=[temp1(:,1),temp1];
temp2=filter(h,1,temp1');
temp2=temp2';
temp3=temp2(:,1:2:2*N);
for c=1:11
error(c,:)=abs(f(i,2:N+1)-temp3(c,:));
#if(key0==3)||(key0==4)
#error(c,:)=abs(f(i,1:N)-temp3(c,:));
#end
end
[min_val,min_idx]=min(error);
for j=1:N
ff(i,j)=temp(min_idx(j),j);
end
end
ff=ff(:,1:N);
fff(1:M,1:2:2*N)=f(1:M,2:N+1);
#if(key0==3)||(key0==4)
#fff(1:M,1:2:2*N)=ff(1:M,1:N);
#end
fff(1:M,2:2:2*N)=ff(1:M,1:N);
f1w=fff';
X=imread('/home/yhy/图片/zebra.bmp');
%f1=X;
%f1=double(f1(:,:,1))/255;
%[M1,M2]=size(f1);
%fff=(fff>=0).*fff;
%error=fff'(1:M1,:,1)-f1;
%MSE=sum(sum(error(3:M1-2,3:M2-2).^2))*255^2/(M1-4)*(M2-4)
%subplot(3,1,3);
imshow(fff');
title(size(fff'));
save data1w f1w;
二、C++代码
代码中增加了一些并行化处理函数，为了加速运行，但是仍然有很大提速空间，可以考虑参照opencv的图像插值源码加速策略去继续做加速优化。这个代码就是直接针对RGB去做处理了。
代码中集成了几种基于opencv的插值方法，自己写了一个自适应插值法，对比了效果（原理参照另一篇博客）。
#include <iostream>
#include <opencv2/opencv.hpp>
#include <vector>
#include <thread>
#include <cmath>
#include <omp.h>
#define PATH "/home/yhy/文档/CPP/Super resolution/SuperResolution/DJI_0044_3.jpg"
#define KEY0 1 //Enter 1 for adaptive cubic,enter 2 for normal cubic,enter for normal lanczos:
using namespace std;
using namespace cv;
Mat imgR;
Mat imgG;
Mat imgB;
int main()
{
Mat treatimg(Mat img,int M,int N);//对一个纬度插值
Mat imgGrey(Mat img,int M,int N);//针对某个channel使用treat函数完成横向和纵向扩充，生成一个扩充后的灰度图像
void imgRED(Mat img,int M,int N);//封装后的GREY,因为多线程thread只能接受void型函数，因此使用全剧变量进行了封装
void imgGREEN(Mat img,int M,int N);
void imgBLUE(Mat img,int M,int N);
vector<Mat> img;
Mat imgrgb = imread(PATH,IMREAD_ANYCOLOR); //读取同目录下的1.jpg文件
int C=imgrgb.channels();
int M=imgrgb.rows;
int N=imgrgb.cols;
int key;
cout<<"image channel is "<<C<<"\n"<<"image dims is "<<M<<" * "<<N<<"\n";//通道数
key=KEY0;
if (imgrgb.channels()==3 and key==1)
{
split(imgrgb,img); //分割
cout<<img[0].channels()<<" is ok"<<"\n";
}
else if(key==1) //灰度图像
{
Mat img[]={imgrgb};
cout<<imgrgb.channels()<<" is ok"<<"\n";
}
if(key==1)
{
//img[0]=imgGrey(img[0],M,N);//R
//img[1]=imgGrey(img[1],M,N);//G
//img[2]=imgGrey(img[2],M,N);//B
thread R(imgRED,img[0],M,N);
thread G(imgGREEN,img[1],M,N);
thread B(imgBLUE,img[2],M,N);
R.join();
G.join();
B.join();
img[0]=imgR;
img[1]=imgG;
img[2]=imgB;
}
if(key==2)
{
resize(imgrgb,imgrgb,Size(),2,2,INTER_CUBIC);
}
if(key==3)
{
resize(imgrgb,imgrgb,Size(),2,2,INTER_LANCZOS4);
}
if(key==1)
{
merge(img,imgrgb); //三通道合成
imgrgb=imgrgb.mul(255);
}
//cout<<imgrgb;
imwrite("2.jpg", imgrgb);
cout<<"treat is ok"<<"image dims is "<<2*M<<" * "<<2*N<<"\n"<<" ";
return 0;
}
Mat treatimg(Mat img,int M,int N)
{
cout<<"Start treat:";
//fliter
float a=-0.5;
int L=256;
float k=1.5;
int l1 = 1,l2 = 1,l3 = 1,l4 = 1;
float error;
float minval[1920];
int minidx[1920]={0};
Mat s = (Mat_<float>(1, 11) << 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1);
Mat h = (Mat_<float>(1, 2) << 0.5, 0.5);
Mat ff =Mat(Size(N,M),CV_32FC(1)) ;ff=0; //M*N的零矩阵
Mat fff =Mat(Size(N*2,M),CV_32FC(1)) ;fff=0; //M*N的零矩阵
Mat temp = Mat(Size(N,11),CV_32FC(1)) ;
Mat temp1 = Mat(Size(2*N+1,11),CV_32FC(1)) ;temp1=0;
Mat temp3 = Mat(Size(N,11),CV_32FC(1)) ;temp3=0;
//Mat error = Mat(Size(N,11),CV_32FC(1)) ;
Mat x=img.colRange(N-2,N).clone(); //沿y轴截取 row：沿x轴截取
flip(x,x,-1); //-1：x与y轴对称变换，0：只x 1：只y
Mat y=img.colRange(0,1).clone(); //不加clone会对元矩阵直接操作
hconcat(y,img,img);
hconcat(img,x,img); //矩阵合并img=[y,img,x]
//cout<<l1*img.at<float>(i,j-1)*(a*pow(s.at<float>(c),3)-2*a*pow(s.at<float>(c),2)+a*s.at<float>(c));
#pragma omp parallel for
for (int i=0;i<=M-1;i++)
{
//#pragma omp parallel for collapse(2)
for (int c=0;c<=10;c++)
//#pragma omp parallel for
for (int j=1;j<=N;j++)
{
temp.at<float>(c,j-1)=l1*img.at<float>(i,j-1)*(a*pow(s.at<float>(c),3)-2*a*pow(s.at<float>(c),2)+a*s.at<float>(c))+
l2*img.at<float>(i,j)*((a+2)*pow(s.at<float>(c),3)-(3+a)*pow(s.at<float>(c),2)+1)+
l3*img.at<float>(i,j+1)*(-(a+2)*pow(s.at<float>(c),3)+(2*a+3)*pow(s.at<float>(c),2)-a*s.at<float>(c))+
l4*img.at<float>(i,j+2)*(-a*pow(s.at<float>(c),3)+a*pow(s.at<float>(c),2));
}
//#pragma omp parallel for collapse(2)
for (int c=0;c<=10;c++)
//#pragma omp parallel for
for (int n=0;n<=N-1;n++)
{
temp1.at<float>(c,2*n)=img.at<float>(i,n+1);
}
//#pragma omp parallel for collapse(2)
for (int c=0;c<=10;c++)
//#pragma omp parallel for
for (int n=0;n<=N-1;n++)
{
temp1.at<float>(c,2*n+1)=temp.at<float>(c,n);
//filter .T
Mat temp2=temp1.colRange(0,1).clone();
//cout<<temp1;
hconcat(temp2,temp1,temp2);
//filter2D(temp2,temp2,-1, 1, Point(-1, -1), 0, 0);
temp3.at<float>(c,n)=temp2.at<float>(c,2*n);
}
for (int n=0;n<=N-1;n++)
{
minval[n]=1;
#pragma omp parallel for
for (int c=0;c<=9;c++)
{
error=abs(img.at<float>(i,n+1)-temp3.at<float>(c,n));
//cout<<error<<" ";
if(minval[n]>=error)
{
minval[n]=error;
minidx[n]=c;
}
}
ff.at<float>(i,n)=temp.at<float>(minidx[n],n);
}
}
ff=ff.colRange(0,N).clone();
//#pragma omp parallel for
#pragma omp parallel for collapse(2)
for(int i=0;i<=M-1;i++)
for(int n=0;n<=N-1;n++)
{
float zero=0;
fff.at<float>(i,2*n)=max(img.at<float>(i,n+1),zero);
fff.at<float>(i,2*n+1)=max(ff.at<float>(i,n),zero); //像素没有负数，计算得到的小负数默认灰度值=0
}
return fff;
}
Mat imgGrey(Mat img,int M,int N)
{
Mat imgcow;
img.convertTo(img, CV_32FC(1), 1 / 255.0);//归一化
imgcow=treatimg(img,M,N); //由M*N变为M*2N
transpose(imgcow,img);
imgcow=treatimg(img,2*N,M);//由M*2N变为2M*2N
transpose(imgcow,imgcow);//图像转置
return imgcow;
}
void imgRED(Mat img,int M,int N)
{
imgR=imgGrey(img,M,N);
}
void imgGREEN(Mat img,int M,int N)
{
imgG=imgGrey(img,M,N);
}
void imgBLUE(Mat img,int M,int N)
{
imgB=imgGrey(img,M,N);
}
三、运行结果
图片一：
原始：
双三次：
lanczos:：
自适应：
图片二：
原始：
双三次：
lanczos:：
自适应：
四、MATLAB总结
同列不同行的两个矩阵，我们可以通过c=[a;b]，进行拼接矩阵，分号是纵向拼接矩阵，
对于同行数不同列数的两个矩阵，可以通过c=[a,d]，进行拼接矩阵，其中逗号的意思是横向拼接矩阵
imread(’/home/yhy/图片/zebra.bmp’) 读取图片
imwrite(E, ‘butterfly_GTME.bmp’) 写图片
imshow§ 显示图片
五、C++总结
推荐C++opencv入门教程：https://blog.csdn.net/qq_42885747/article/details/109729532
const传入不变量，如果传入变量前加&表示送入地址，调用时不用增加额外操作，少用内存
图像处理中变量命名
imageTemp
imageInput
mageOutput
imageCol
imageRow
Mat y=img.colRange(0,1).clone();
连续几行的复制
需要注意的是，上面两个函数返回的子矩阵是对原矩阵的引用。即如果修改子矩阵r_range，原矩阵相应地也会被改变。如果不想改变原矩阵中的值，即想要获取一个拷贝的子矩阵，可以使用clone函数：
扫一扫
839
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
