鱼眼拼接之SIFT提取特征点
2015-06-10 17:31
429 查看
之前写的那两篇 是我没有完全理解SIFT理论部分 后来看了这位女侠的http://blog.csdn.net/abcjennifer/article/details/7639681终于弄明白了
%这部分是之前写过的鱼眼提取有效区域及经度坐标校正
A=imread('F:\orl_zhifangtu\s1.jpg');
T=40;
[A,R]=kuaisusaomiao(A,T);
[m,n,H]=size(A);
C=[];
x=n/2;
y=m/2;
for u=1:m
for v=1:n
i=u;
j=round(sqrt(R^2-(y-u)^2)*(v-x)/R+x);
if(R^2-(y-u)^2<0)
continue;
end
C(u,v,1)=A(i,j,1);
C(u,v,2)=A(i,j,2);
C(u,v,3)=A(i,j,3);
end
end
C=uint8(C);
这是原图
这是校正然后灰度化以后的图
接下来是SIFT高斯卷积 这里之前我用的http://blog.csdn.net/v_JULY_v/article/details/6245939这个理论部分 按照他的高斯函数sigma第一层第一幅图取值是1.6 第二幅图是k*1.6 第三幅是K^2*1.6 第四幅是K^3*1.6 第五幅是K^4*1.6（第一层五幅图的高斯函数已经完毕 接下来第二层）;第二层第一幅图是2*1.6 第二幅图是2*K^2*1.6（是上一层第二幅的两倍） 第三幅是2*K^3*1.6 以此类推。。。第三层第一幅图是2*2*1.6 第二幅图是2*2*K*1.6（即是上一层的两倍）以此类推。。。。但完完全全按照这个大神编
后面会有问题 就是比较极值显示时 会一个点都没有 （我忘记了当时是一个点都没有还是只有肉眼可见的零星的几个点 反正不是我们要的特征点） 所以后来看这个女侠解说的理论 尺度连续性那里总算看懂了 但是编的时候如果完全按照那个参数sigma编 还是会有一样的问题
7，0.7
所以我总共编了三层金字塔，每一层有六幅图像，每一次的sigma不是按照k变化的 而是随意取的0.3，0.4，0.5，0.6，0.7，0.8这六个参数 可得到六个高斯函数 然后对原图进行卷积 得到 第一层高斯卷积后的图
%先灰度化
C=rgb2gray(C);
C=double(C);
[m1,n1,h1]=size(C);
k=2^(1/3);
threshold=3;
h11=fspecial('gaussian',[5 5],0.3);
C11=imfilter(C,h11,'conv');
h12=fspecial('gaussian',[5 5],0.4);
C12=imfilter(C,h12,'conv');
h13=fspecial('gaussian',[5 5],0.5);
C13=imfilter(C,h13,'conv');
h14=fspecial('gaussian',[5 5],0.6);
C14=imfilter(C,h14,'conv');
h15=fspecial('gaussian',[5 5],0.7);
C15=imfilter(C,h15,'conv');
h16=fspecial('gaussian',[5 5],0.8);
C16=imfilter(C,h16,'conv');
%以上就得到了第一层的六幅高斯卷积后的图C11，C12，C13,C14，C15，C16
然后两两相减 得到第一层的DOG差分空间后的图 只有五幅 D11，D12，D13，D14，D15
D11=C11-C12;
D12=C12-C13;
D13=C13-C14;
D14=C14-C15;
D15=C15-C16;
%接下来就是对第一层进行极值比较 初步得到候选的特征点 只有D12 D13 D14前后都有图 所以才可进行比较 首尾即D11，D15没办法得到26个特征点 就不能进行比较 所以我把D12 D13 D14这三幅图得到的特征点 分别装在E11 E12 E13里
E11=zeros(m1,n1);
for i=3:m1-2
for j=3:n1-2
if(D12(i,j)>D12(i-1,j-1) && D12(i,j)>D12(i,j-1) && D12(i,j)>D12(i+1,j-1) && D12(i,j)>D12(i-1,j) && D12(i,j)>D12(i+1,j) && D12(i,j)>D12(i-1,j+1) && D12(i,j)>D12(i,j+1) && D12(i,j)>D12(i+1,j+1))
if(D12(i,j)>D11(i,j) && D12(i,j)>D11(i-1,j-1) && D12(i,j)>D11(i,j-1) && D12(i,j)>D11(i+1,j-1) && D12(i,j)>D11(i-1,j) && D12(i,j)>D11(i+1,j) && D12(i,j)>D11(i-1,j+1) && D12(i,j)>D11(i,j+1) && D12(i,j)>D11(i+1,j+1))
if(D12(i,j)>D13(i,j) && D12(i,j)>D13(i-1,j-1) && D12(i,j)>D13(i,j-1) && D12(i,j)>D13(i+1,j-1) && D12(i,j)>D13(i-1,j) && D12(i,j)>D13(i+1,j) && D12(i,j)>D13(i-1,j+1) && D12(i,j)>D13(i,j+1) && D12(i,j)>D13(i+1,j+1))
if(D12(i,j)>threshold)
E11(i,j)=1;
end
end
end
elseif(D12(i,j)<D12(i-1,j-1) && D12(i,j)<D12(i,j-1) && D12(i,j)<D12(i+1,j-1) && D12(i,j)<D12(i-1,j) && D12(i,j)<D12(i+1,j) && D12(i,j)<D12(i-1,j+1) && D12(i,j)<D12(i,j+1) && D12(i,j)<D12(i+1,j+1))
if(D12(i,j)<D11(i,j) && D12(i,j)<D11(i-1,j-1) && D12(i,j)<D11(i,j-1) && D12(i,j)<D11(i+1,j-1) && D12(i,j)<D11(i-1,j) && D12(i,j)<D11(i+1,j) && D12(i,j)<D11(i-1,j+1) && D12(i,j)<D11(i,j+1) && D12(i,j)<D11(i+1,j+1))
if(D12(i,j)<D13(i,j) && D12(i,j)<D13(i-1,j-1) && D12(i,j)<D13(i,j-1) && D12(i,j)<D13(i+1,j-1) && D12(i,j)<D13(i-1,j) && D12(i,j)<D13(i+1,j) && D12(i,j)<D13(i-1,j+1) && D12(i,j)<D13(i,j+1) && D12(i,j)<D13(i+1,j+1))
if(D12(i,j)<-threshold)
E11(i,j)=-1;
end
end
end
end
end
end
E12=zeros(m1,n1);
for i=2:m1-1
for j=2:n1-1
if(D13(i,j)>D13(i-1,j-1) && D13(i,j)>D13(i,j-1) && D13(i,j)>D13(i+1,j-1) && D13(i,j)>D13(i-1,j) && D13(i,j)>D13(i+1,j) && D13(i,j)>D13(i-1,j+1) && D13(i,j)>D12(i,j+1) && D13(i,j)>D13(i+1,j+1))
if(D13(i,j)>D12(i,j) && D13(i,j)>D12(i-1,j-1) && D13(i,j)>D12(i,j-1) && D13(i,j)>D12(i+1,j-1) && D13(i,j)>D12(i-1,j) && D13(i,j)>D12(i+1,j) && D13(i,j)>D12(i-1,j+1) && D13(i,j)>D12(i,j+1) && D13(i,j)>D12(i+1,j+1))
if(D13(i,j)>D14(i,j) && D13(i,j)>D14(i-1,j-1) && D13(i,j)>D14(i,j-1) && D13(i,j)>D14(i+1,j-1) && D13(i,j)>D14(i-1,j) && D13(i,j)>D14(i+1,j) && D13(i,j)>D14(i-1,j+1) && D13(i,j)>D14(i,j+1) && D13(i,j)>D14(i+1,j+1))
if(D13(i,j)>threshold)
E12(i,j)=1;
end
end
end
elseif(D13(i,j)<D13(i-1,j-1) && D13(i,j)<D13(i,j-1) && D13(i,j)<D13(i+1,j-1) && D13(i,j)<D13(i-1,j) && D13(i,j)<D13(i+1,j) && D13(i,j)<D13(i-1,j+1) && D13(i,j)<D12(i,j+1) && D13(i,j)<D13(i+1,j+1))
if(D13(i,j)<D12(i,j) && D13(i,j)<D12(i-1,j-1) && D13(i,j)<D12(i,j-1) && D13(i,j)<D12(i+1,j-1) && D13(i,j)<D12(i-1,j) && D13(i,j)<D12(i+1,j) && D13(i,j)<D12(i-1,j+1) && D13(i,j)<D12(i,j+1) && D13(i,j)<D12(i+1,j+1))
if(D13(i,j)<D14(i,j) && D13(i,j)<D14(i-1,j-1) && D13(i,j)<D14(i,j-1) && D13(i,j)<D14(i+1,j-1) && D13(i,j)<D14(i-1,j) && D13(i,j)<D14(i+1,j) && D13(i,j)<D14(i-1,j+1) && D13(i,j)<D14(i,j+1) && D13(i,j)<D14(i+1,j+1))
if(D13(i,j)<-threshold)
E12(i,j)=-1;
end
end
end
end
end
end
E13=zeros(m1,n1);
for i=2:m1-1
for j=2:n1-1
if(D14(i,j)>D14(i-1,j-1) && D14(i,j)>D14(i,j-1) && D14(i,j)>D14(i+1,j-1) && D14(i,j)>D14(i-1,j) && D14(i,j)>D14(i+1,j) && D14(i,j)>D14(i-1,j+1) && D14(i,j)>D14(i,j+1) && D14(i,j)>D14(i+1,j+1))
if(D14(i,j)>D13(i,j) && D14(i,j)>D13(i-1,j-1) && D14(i,j)>D13(i,j-1) && D14(i,j)>D13(i+1,j-1) && D14(i,j)>D13(i-1,j) && D14(i,j)>D13(i+1,j) && D14(i,j)>D13(i-1,j+1) && D14(i,j)>D13(i,j+1) && D14(i,j)>D13(i+1,j+1))
if(D14(i,j)>D15(i,j) && D14(i,j)>D15(i-1,j-1) && D14(i,j)>D15(i,j-1) && D14(i,j)>D15(i+1,j-1) && D14(i,j)>D15(i-1,j) && D14(i,j)>D15(i+1,j) && D14(i,j)>D15(i-1,j+1) && D14(i,j)>D15(i,j+1) && D14(i,j)>D15(i+1,j+1))
if(D14(i,j)>threshold)
E13(i,j)=1;
end
end
end
elseif(D14(i,j)<D14(i-1,j-1) && D14(i,j)<D14(i,j-1) && D14(i,j)<D14(i+1,j-1) && D14(i,j)<D14(i-1,j) && D14(i,j)<D14(i+1,j) && D14(i,j)<D14(i-1,j+1) && D14(i,j)<D14(i,j+1) && D14(i,j)<D14(i+1,j+1))
if(D14(i,j)<D13(i,j) && D14(i,j)<D13(i-1,j-1) && D14(i,j)<D13(i,j-1) && D14(i,j)<D13(i+1,j-1) && D14(i,j)<D13(i-1,j) && D14(i,j)<D13(i+1,j) && D14(i,j)<D13(i-1,j+1) && D14(i,j)<D13(i,j+1) && D14(i,j)<D13(i+1,j+1))
if(D14(i,j)<D15(i,j) && D14(i,j)<D15(i-1,j-1) && D14(i,j)<D15(i,j-1) && D14(i,j)<D15(i+1,j-1) && D14(i,j)<D15(i-1,j) && D14(i,j)<D15(i+1,j) && D14(i,j)<D15(i-1,j+1) && D14(i,j)<D15(i,j+1) && D14(i,j)<D15(i+1,j+1))
if(D14(i,j)>threshold)
E13(i,j)=-1;
end
end
end
end
end
end
这样就得到了初步得到了第一层的候选特征点
figure
>> imshow(D11,[0 1])
>> figure
>> imshow(D12,[0 1])
>> figure
>> imshow(D13,[0 1])
>> figure
>> imshow(D14,[0 1])
>> figure
>> imshow(D15,[0 1])
这是D11
这是D12
这是D13
这是D14
这是D15
%下面显示E11 即极大极小值 特征点
figure
>> imshow(E11,[-1 1])
>> figure
>> imshow(E12,[-1 1])
>> figure
>> imshow(E13,[-1 1])
白色表示极大值 黑色表示极小值
但是E12 E13这两幅 出来就看不到特征点了 我在想 是不是 一层都只有一幅有特征点的图出来 后面的两幅图都看不到特征点
下面是E12 和 E13
可以看到
两幅都没有特征点出来 好奇怪明明都是26个像素比较 容我想想 不过既然只有E11有特征点出来 我们就暂时第一层就用E11吧 至少出来了一个
第二层 第三层跟上面步骤道理一样：
m2=round(m1/2);
n2=round(n1/2);
C20=imresize(C11,[m2 n2]);
h21=fspecial('gaussian',[5 5],k*0.3);
C21=imfilter(C20,h21,'conv');
h22=fspecial('gaussian',[5 5],k*0.4);
C22=imfilter(C20,h22,'conv');
h23=fspecial('gaussian',[5 5],k*0.5);
C23=imfilter(C20,h23,'conv');
h24=fspecial('gaussian',[5 5],k*0.6);
C24=imfilter(C20,h24,'conv');
h25=fspecial('gaussian',[5 5],k*0.7);
C25=imfilter(C20,h25,'conv');
h26=fspecial('gaussian',[5 5],k*0.8);
C26=imfilter(C20,h26,'conv');
m3=round(m2/2);
n3=round(n2/2);
C30=imresize(C21,[m3 n3]);
h31=fspecial('gaussian',[5 5],k^2*0.3);
C31=imfilter(C30,h31,'conv');
h32=fspecial('gaussian',[5 5],k^2*0.4);
C32=imfilter(C30,h32,'conv');
h33=fspecial('gaussian',[5 5],k^2*0.5);
C33=imfilter(C30,h33,'conv');
h34=fspecial('gaussian',[5 5],k^2*0.6);
C34=imfilter(C30,h34,'conv');
h35=fspecial('gaussian',[5 5],k^2*0.7);
C35=imfilter(C30,h35,'conv');
h36=fspecial('gaussian',[5 5],k^2*0.8);
C36=imfilter(C30,h36,'conv');
>> D21=C21-C22;
D22=C22-C23;
D23=C23-C24;
D24=C24-C25;
D25=C25-C26;
D31=C31-C32;
D32=C32-C33;
D33=C33-C34;
D34=C34-C35;
D35=C35-C36;
>> figure
>> imshow(D21,[0 1])
这是第二层的第一幅DOG高斯尺度空间的图 依次可以显示第二层第三层的所有图。。。
比较第二层的极值点
E21=zeros(m2,n2);
for i=2:m2-1
for j=2:n2-1
if(D22(i,j)>D22(i-1,j-1) && D22(i,j)>D22(i,j-1) && D22(i,j)>D22(i+1,j-1) && D22(i,j)>D22(i-1,j) && D22(i,j)>D22(i+1,j) && D22(i,j)>D22(i-1,j+1) && D22(i,j)>D22(i,j+1) && D22(i,j)>D22(i+1,j+1))
if(D22(i,j)>D21(i,j) && D22(i,j)>D21(i-1,j-1) && D22(i,j)>D21(i,j-1) && D22(i,j)>D21(i+1,j-1) && D22(i,j)>D21(i-1,j) && D22(i,j)>D21(i+1,j) && D22(i,j)>D21(i-1,j+1) && D22(i,j)>D21(i,j+1) && D22(i,j)>D21(i+1,j+1))
if(D22(i,j)>D23(i,j) && D22(i,j)>D23(i-1,j-1) && D22(i,j)>D23(i,j-1) && D22(i,j)>D23(i+1,j-1) && D22(i,j)>D23(i-1,j) && D22(i,j)>D23(i+1,j) && D22(i,j)>D23(i-1,j+1) && D22(i,j)>D23(i,j+1) && D22(i,j)>D23(i+1,j+1))
if(D22(i,j)>threshold)
E21(i,j)=1;
end
end
end
elseif(D22(i,j)<D22(i-1,j-1) && D22(i,j)<D22(i,j-1) && D22(i,j)<D22(i+1,j-1) && D22(i,j)<D22(i-1,j) && D22(i,j)<D22(i+1,j) && D22(i,j)<D22(i-1,j+1) && D22(i,j)<D22(i,j+1) && D22(i,j)<D22(i+1,j+1))
if(D22(i,j)<D21(i,j) && D22(i,j)<D21(i-1,j-1) && D22(i,j)<D21(i,j-1) && D22(i,j)<D21(i+1,j-1) && D22(i,j)<D21(i-1,j) && D22(i,j)<D21(i+1,j) && D22(i,j)<D21(i-1,j+1) && D22(i,j)<D21(i,j+1) && D22(i,j)<D21(i+1,j+1))
if(D22(i,j)<D23(i,j) && D22(i,j)<D23(i-1,j-1) && D22(i,j)<D23(i,j-1) && D22(i,j)<D23(i+1,j-1) && D22(i,j)<D23(i-1,j) && D22(i,j)<D23(i+1,j) && D22(i,j)<D23(i-1,j+1) && D22(i,j)<D23(i,j+1) && D22(i,j)<D23(i+1,j+1))
if(D22(i,j)<-threshold)
E21(i,j)=-1;
end
end
end
end
end
end
E22=zeros(m2,n2);
for i=2:m2-1
for j=2:n2-1
if(D23(i,j)>D23(i-1,j-1) && D23(i,j)>D23(i,j-1) && D23(i,j)>D23(i+1,j-1) && D23(i,j)>D23(i-1,j) && D23(i,j)>D23(i+1,j) && D23(i,j)>D23(i-1,j+1) && D23(i,j)>D23(i,j+1) && D23(i,j)>D23(i+1,j+1))
if(D23(i,j)>D22(i,j) && D23(i,j)>D22(i-1,j-1) && D23(i,j)>D22(i,j-1) && D23(i,j)>D22(i+1,j-1) && D23(i,j)>D22(i-1,j) && D23(i,j)>D22(i+1,j) && D23(i,j)>D22(i-1,j+1) && D23(i,j)>D22(i,j+1) && D23(i,j)>D22(i+1,j+1))
if(D23(i,j)>D24(i,j) && D23(i,j)>D24(i-1,j-1) && D23(i,j)>D24(i,j-1) && D23(i,j)>D24(i+1,j-1) && D23(i,j)>D24(i-1,j) && D23(i,j)>D24(i+1,j) && D23(i,j)>D24(i-1,j+1) && D23(i,j)>D24(i,j+1) && D23(i,j)>D24(i+1,j+1))
if(D23(i,j)>threshold)
E22(i,j)=1;
end
end
end
elseif(D23(i,j)<D23(i-1,j-1) && D23(i,j)<D23(i,j-1) && D23(i,j)<D23(i+1,j-1) && D23(i,j)<D23(i-1,j) && D23(i,j)<D23(i+1,j) && D23(i,j)<D23(i-1,j+1) && D23(i,j)<D23(i,j+1) && D23(i,j)<D23(i+1,j+1))
if(D23(i,j)<D22(i,j) && D23(i,j)<D22(i-1,j-1) && D23(i,j)<D22(i,j-1) && D23(i,j)<D22(i+1,j-1) && D23(i,j)<D22(i-1,j) && D23(i,j)<D22(i+1,j) && D23(i,j)<D22(i-1,j+1) && D23(i,j)<D22(i,j+1) && D23(i,j)<D22(i+1,j+1))
if(D23(i,j)<D24(i,j) && D23(i,j)<D24(i-1,j-1) && D23(i,j)<D24(i,j-1) && D23(i,j)<D24(i+1,j-1) && D23(i,j)<D24(i-1,j) && D23(i,j)<D24(i+1,j) && D23(i,j)<D24(i-1,j+1) && D23(i,j)<D24(i,j+1) && D23(i,j)<D24(i+1,j+1))
if(D23(i,j)<-threshold)
E22(i,j)=-1;
end
end
end
end
end
end
E23=zeros(m2,n2);
for i=2:m2-1
for j=2:n2-1
if(D24(i,j)>D24(i-1,j-1) && D24(i,j)>D24(i,j-1) && D24(i,j)>D24(i+1,j-1) && D24(i,j)>D24(i-1,j) && D24(i,j)>D24(i+1,j) && D24(i,j)>D24(i-1,j+1) && D24(i,j)>D24(i,j+1) && D24(i,j)>D24(i+1,j+1))
if(D24(i,j)>D23(i,j) && D24(i,j)>D23(i-1,j-1) && D24(i,j)>D23(i,j-1) && D24(i,j)>D23(i+1,j-1) && D24(i,j)>D23(i-1,j) && D24(i,j)>D23(i+1,j) && D24(i,j)>D23(i-1,j+1) && D24(i,j)>D23(i,j+1) && D24(i,j)>D23(i+1,j+1))
if(D24(i,j)>D25(i,j) && D24(i,j)>D25(i-1,j-1) && D24(i,j)>D25(i,j-1) && D24(i,j)>D25(i+1,j-1) && D24(i,j)>D25(i-1,j) && D24(i,j)>D25(i+1,j) && D24(i,j)>D25(i-1,j+1) && D24(i,j)>D25(i,j+1) && D24(i,j)>D25(i+1,j+1))
if(D24(i,j)>threshold)
E23(i,j)=1;
end
end
end
elseif(D24(i,j)<D24(i-1,j-1) && D24(i,j)<D24(i,j-1) && D24(i,j)<D24(i+1,j-1) && D24(i,j)<D24(i-1,j) && D24(i,j)<D24(i+1,j) && D24(i,j)<D24(i-1,j+1) && D24(i,j)<D24(i,j+1) && D24(i,j)<D24(i+1,j+1))
if(D24(i,j)<D23(i,j) && D24(i,j)<D23(i-1,j-1) && D24(i,j)<D23(i,j-1) && D24(i,j)<D23(i+1,j-1) && D24(i,j)<D23(i-1,j) && D24(i,j)<D23(i+1,j) && D24(i,j)<D23(i-1,j+1) && D24(i,j)<D23(i,j+1) && D24(i,j)<D23(i+1,j+1))
if(D24(i,j)<D25(i,j) && D24(i,j)<D25(i-1,j-1) && D24(i,j)<D25(i,j-1) && D24(i,j)<D25(i+1,j-1) && D24(i,j)<D25(i-1,j) && D24(i,j)<D25(i+1,j) && D24(i,j)<D25(i-1,j+1) && D24(i,j)<D25(i,j+1) && D24(i,j)<D25(i+1,j+1))
if(D24(i,j)<-threshold)
E23(i,j)=-1;
end
end
end
end
end
end
下面显示E21
figure
>> imshow(E21,[-1 1])
不知你们看到这幅图的四五个黑白点没 就是极值点 我自己也是醉了
后面第三层的就不显示了 我估计就是一片灰色 连这四五个点都不会给我显示出来了。。。
经过上面的分析以及显示的图像 我们知道其实只需要第一层就可以得到候选特征点 所以我只要E11 完整的程序如下所示：
A=imread('F:\orl_zhifangtu\s1.jpg');
T=40;
[A,R]=kuaisusaomiao(A,T);
[m,n,H]=size(A);
C=[];
x=n/2;
y=m/2;
for u=1:m
for v=1:n
i=u;
j=round(sqrt(R^2-(y-u)^2)*(v-x)/R+x);
if(R^2-(y-u)^2<0)
continue;
end
C(u,v,1)=A(i,j,1);
C(u,v,2)=A(i,j,2);
C(u,v,3)=A(i,j,3);
end
end
C=uint8(C);
C=rgb2gray(C);
H=C;
C=double(C);
[m1,n1,h1]=size(C);
k=2^(1/3);
threshold=3;
h11=fspecial('gaussian',[5 5],0.3);
C11=imfilter(C,h11,'conv');
h12=fspecial('gaussian',[5 5],0.4);
C12=imfilter(C,h12,'conv');
h13=fspecial('gaussian',[5 5],0.5);
C13=imfilter(C,h13,'conv');
h14=fspecial('gaussian',[5 5],0.6);
C14=imfilter(C,h14,'conv');
h15=fspecial('gaussian',[5 5],0.7);
C15=imfilter(C,h15,'conv');
h16=fspecial('gaussian',[5 5],0.8);
C16=imfilter(C,h16,'conv');
D11=C11-C12;
D12=C12-C13;
D13=C13-C14;
D14=C14-C15;
D15=C15-C16;
E11=zeros(m1,n1);
for i=3:m1-2
for j=3:n1-2
if(D12(i,j)>D12(i-1,j-1) && D12(i,j)>D12(i,j-1) && D12(i,j)>D12(i+1,j-1) && D12(i,j)>D12(i-1,j) && D12(i,j)>D12(i+1,j) && D12(i,j)>D12(i-1,j+1) && D12(i,j)>D12(i,j+1) && D12(i,j)>D12(i+1,j+1))
if(D12(i,j)>D11(i,j) && D12(i,j)>D11(i-1,j-1) && D12(i,j)>D11(i,j-1) && D12(i,j)>D11(i+1,j-1) && D12(i,j)>D11(i-1,j) && D12(i,j)>D11(i+1,j) && D12(i,j)>D11(i-1,j+1) && D12(i,j)>D11(i,j+1) && D12(i,j)>D11(i+1,j+1))
if(D12(i,j)>D13(i,j) && D12(i,j)>D13(i-1,j-1) && D12(i,j)>D13(i,j-1) && D12(i,j)>D13(i+1,j-1) && D12(i,j)>D13(i-1,j) && D12(i,j)>D13(i+1,j) && D12(i,j)>D13(i-1,j+1) && D12(i,j)>D13(i,j+1) && D12(i,j)>D13(i+1,j+1))
if(D12(i,j)>threshold)
E11(i,j)=1;
end
end
end
elseif(D12(i,j)<D12(i-1,j-1) && D12(i,j)<D12(i,j-1) && D12(i,j)<D12(i+1,j-1) && D12(i,j)<D12(i-1,j) && D12(i,j)<D12(i+1,j) && D12(i,j)<D12(i-1,j+1) && D12(i,j)<D12(i,j+1) && D12(i,j)<D12(i+1,j+1))
if(D12(i,j)<D11(i,j) && D12(i,j)<D11(i-1,j-1) && D12(i,j)<D11(i,j-1) && D12(i,j)<D11(i+1,j-1) && D12(i,j)<D11(i-1,j) && D12(i,j)<D11(i+1,j) && D12(i,j)<D11(i-1,j+1) && D12(i,j)<D11(i,j+1) && D12(i,j)<D11(i+1,j+1))
if(D12(i,j)<D13(i,j) && D12(i,j)<D13(i-1,j-1) && D12(i,j)<D13(i,j-1) && D12(i,j)<D13(i+1,j-1) && D12(i,j)<D13(i-1,j) && D12(i,j)<D13(i+1,j) && D12(i,j)<D13(i-1,j+1) && D12(i,j)<D13(i,j+1) && D12(i,j)<D13(i+1,j+1))
if(D12(i,j)<-threshold)
E11(i,j)=-1;
end
end
end
end
end
end
%下面这部分是把特征点显示在一个原图中 而不是像E11一样显示在一个空矩阵中
[X11,Y11]=find(E11==1 | E11==-1);
imshow(H)
hold on;
plot(Y11,X11,'rx');（这里敲键盘的时候千万别敲错啊 我之前敲错了 换了X11 Y11的位置 结果出来是个亮瞎我铝合金狗眼的图 我自己想了两天 还以为是我哪里编错了）
这就是以上SIFT程序得到的特征点结果 勉强可以吧
后来搜了http://blog.csdn.net/abcjennifer/article/details/7639681她写的MATLAB程序 我真的有一巴掌拍死自己的冲动 人家写得多么简洁 而我的。。。让我自己一刀了断。。。所以你们还是直接下她的程序比较好 短小精悍 又容易懂
我写这个就当是自己的学习笔记好了。。。用来以后 自己边吐血边回忆。。。
%这部分是之前写过的鱼眼提取有效区域及经度坐标校正
A=imread('F:\orl_zhifangtu\s1.jpg');
T=40;
[A,R]=kuaisusaomiao(A,T);
[m,n,H]=size(A);
C=[];
x=n/2;
y=m/2;
for u=1:m
for v=1:n
i=u;
j=round(sqrt(R^2-(y-u)^2)*(v-x)/R+x);
if(R^2-(y-u)^2<0)
continue;
end
C(u,v,1)=A(i,j,1);
C(u,v,2)=A(i,j,2);
C(u,v,3)=A(i,j,3);
end
end
C=uint8(C);
这是原图
这是校正然后灰度化以后的图
接下来是SIFT高斯卷积 这里之前我用的http://blog.csdn.net/v_JULY_v/article/details/6245939这个理论部分 按照他的高斯函数sigma第一层第一幅图取值是1.6 第二幅图是k*1.6 第三幅是K^2*1.6 第四幅是K^3*1.6 第五幅是K^4*1.6（第一层五幅图的高斯函数已经完毕 接下来第二层）;第二层第一幅图是2*1.6 第二幅图是2*K^2*1.6（是上一层第二幅的两倍） 第三幅是2*K^3*1.6 以此类推。。。第三层第一幅图是2*2*1.6 第二幅图是2*2*K*1.6（即是上一层的两倍）以此类推。。。。但完完全全按照这个大神编
后面会有问题 就是比较极值显示时 会一个点都没有 （我忘记了当时是一个点都没有还是只有肉眼可见的零星的几个点 反正不是我们要的特征点） 所以后来看这个女侠解说的理论 尺度连续性那里总算看懂了 但是编的时候如果完全按照那个参数sigma编 还是会有一样的问题
7，0.7
所以我总共编了三层金字塔，每一层有六幅图像，每一次的sigma不是按照k变化的 而是随意取的0.3，0.4，0.5，0.6，0.7，0.8这六个参数 可得到六个高斯函数 然后对原图进行卷积 得到 第一层高斯卷积后的图
%先灰度化
C=rgb2gray(C);
C=double(C);
[m1,n1,h1]=size(C);
k=2^(1/3);
threshold=3;
h11=fspecial('gaussian',[5 5],0.3);
C11=imfilter(C,h11,'conv');
h12=fspecial('gaussian',[5 5],0.4);
C12=imfilter(C,h12,'conv');
h13=fspecial('gaussian',[5 5],0.5);
C13=imfilter(C,h13,'conv');
h14=fspecial('gaussian',[5 5],0.6);
C14=imfilter(C,h14,'conv');
h15=fspecial('gaussian',[5 5],0.7);
C15=imfilter(C,h15,'conv');
h16=fspecial('gaussian',[5 5],0.8);
C16=imfilter(C,h16,'conv');
%以上就得到了第一层的六幅高斯卷积后的图C11，C12，C13,C14，C15，C16
然后两两相减 得到第一层的DOG差分空间后的图 只有五幅 D11，D12，D13，D14，D15
D11=C11-C12;
D12=C12-C13;
D13=C13-C14;
D14=C14-C15;
D15=C15-C16;
%接下来就是对第一层进行极值比较 初步得到候选的特征点 只有D12 D13 D14前后都有图 所以才可进行比较 首尾即D11，D15没办法得到26个特征点 就不能进行比较 所以我把D12 D13 D14这三幅图得到的特征点 分别装在E11 E12 E13里
E11=zeros(m1,n1);
for i=3:m1-2
for j=3:n1-2
if(D12(i,j)>D12(i-1,j-1) && D12(i,j)>D12(i,j-1) && D12(i,j)>D12(i+1,j-1) && D12(i,j)>D12(i-1,j) && D12(i,j)>D12(i+1,j) && D12(i,j)>D12(i-1,j+1) && D12(i,j)>D12(i,j+1) && D12(i,j)>D12(i+1,j+1))
if(D12(i,j)>D11(i,j) && D12(i,j)>D11(i-1,j-1) && D12(i,j)>D11(i,j-1) && D12(i,j)>D11(i+1,j-1) && D12(i,j)>D11(i-1,j) && D12(i,j)>D11(i+1,j) && D12(i,j)>D11(i-1,j+1) && D12(i,j)>D11(i,j+1) && D12(i,j)>D11(i+1,j+1))
if(D12(i,j)>D13(i,j) && D12(i,j)>D13(i-1,j-1) && D12(i,j)>D13(i,j-1) && D12(i,j)>D13(i+1,j-1) && D12(i,j)>D13(i-1,j) && D12(i,j)>D13(i+1,j) && D12(i,j)>D13(i-1,j+1) && D12(i,j)>D13(i,j+1) && D12(i,j)>D13(i+1,j+1))
if(D12(i,j)>threshold)
E11(i,j)=1;
end
end
end
elseif(D12(i,j)<D12(i-1,j-1) && D12(i,j)<D12(i,j-1) && D12(i,j)<D12(i+1,j-1) && D12(i,j)<D12(i-1,j) && D12(i,j)<D12(i+1,j) && D12(i,j)<D12(i-1,j+1) && D12(i,j)<D12(i,j+1) && D12(i,j)<D12(i+1,j+1))
if(D12(i,j)<D11(i,j) && D12(i,j)<D11(i-1,j-1) && D12(i,j)<D11(i,j-1) && D12(i,j)<D11(i+1,j-1) && D12(i,j)<D11(i-1,j) && D12(i,j)<D11(i+1,j) && D12(i,j)<D11(i-1,j+1) && D12(i,j)<D11(i,j+1) && D12(i,j)<D11(i+1,j+1))
if(D12(i,j)<D13(i,j) && D12(i,j)<D13(i-1,j-1) && D12(i,j)<D13(i,j-1) && D12(i,j)<D13(i+1,j-1) && D12(i,j)<D13(i-1,j) && D12(i,j)<D13(i+1,j) && D12(i,j)<D13(i-1,j+1) && D12(i,j)<D13(i,j+1) && D12(i,j)<D13(i+1,j+1))
if(D12(i,j)<-threshold)
E11(i,j)=-1;
end
end
end
end
end
end
E12=zeros(m1,n1);
for i=2:m1-1
for j=2:n1-1
if(D13(i,j)>D13(i-1,j-1) && D13(i,j)>D13(i,j-1) && D13(i,j)>D13(i+1,j-1) && D13(i,j)>D13(i-1,j) && D13(i,j)>D13(i+1,j) && D13(i,j)>D13(i-1,j+1) && D13(i,j)>D12(i,j+1) && D13(i,j)>D13(i+1,j+1))
if(D13(i,j)>D12(i,j) && D13(i,j)>D12(i-1,j-1) && D13(i,j)>D12(i,j-1) && D13(i,j)>D12(i+1,j-1) && D13(i,j)>D12(i-1,j) && D13(i,j)>D12(i+1,j) && D13(i,j)>D12(i-1,j+1) && D13(i,j)>D12(i,j+1) && D13(i,j)>D12(i+1,j+1))
if(D13(i,j)>D14(i,j) && D13(i,j)>D14(i-1,j-1) && D13(i,j)>D14(i,j-1) && D13(i,j)>D14(i+1,j-1) && D13(i,j)>D14(i-1,j) && D13(i,j)>D14(i+1,j) && D13(i,j)>D14(i-1,j+1) && D13(i,j)>D14(i,j+1) && D13(i,j)>D14(i+1,j+1))
if(D13(i,j)>threshold)
E12(i,j)=1;
end
end
end
elseif(D13(i,j)<D13(i-1,j-1) && D13(i,j)<D13(i,j-1) && D13(i,j)<D13(i+1,j-1) && D13(i,j)<D13(i-1,j) && D13(i,j)<D13(i+1,j) && D13(i,j)<D13(i-1,j+1) && D13(i,j)<D12(i,j+1) && D13(i,j)<D13(i+1,j+1))
if(D13(i,j)<D12(i,j) && D13(i,j)<D12(i-1,j-1) && D13(i,j)<D12(i,j-1) && D13(i,j)<D12(i+1,j-1) && D13(i,j)<D12(i-1,j) && D13(i,j)<D12(i+1,j) && D13(i,j)<D12(i-1,j+1) && D13(i,j)<D12(i,j+1) && D13(i,j)<D12(i+1,j+1))
if(D13(i,j)<D14(i,j) && D13(i,j)<D14(i-1,j-1) && D13(i,j)<D14(i,j-1) && D13(i,j)<D14(i+1,j-1) && D13(i,j)<D14(i-1,j) && D13(i,j)<D14(i+1,j) && D13(i,j)<D14(i-1,j+1) && D13(i,j)<D14(i,j+1) && D13(i,j)<D14(i+1,j+1))
if(D13(i,j)<-threshold)
E12(i,j)=-1;
end
end
end
end
end
end
E13=zeros(m1,n1);
for i=2:m1-1
for j=2:n1-1
if(D14(i,j)>D14(i-1,j-1) && D14(i,j)>D14(i,j-1) && D14(i,j)>D14(i+1,j-1) && D14(i,j)>D14(i-1,j) && D14(i,j)>D14(i+1,j) && D14(i,j)>D14(i-1,j+1) && D14(i,j)>D14(i,j+1) && D14(i,j)>D14(i+1,j+1))
if(D14(i,j)>D13(i,j) && D14(i,j)>D13(i-1,j-1) && D14(i,j)>D13(i,j-1) && D14(i,j)>D13(i+1,j-1) && D14(i,j)>D13(i-1,j) && D14(i,j)>D13(i+1,j) && D14(i,j)>D13(i-1,j+1) && D14(i,j)>D13(i,j+1) && D14(i,j)>D13(i+1,j+1))
if(D14(i,j)>D15(i,j) && D14(i,j)>D15(i-1,j-1) && D14(i,j)>D15(i,j-1) && D14(i,j)>D15(i+1,j-1) && D14(i,j)>D15(i-1,j) && D14(i,j)>D15(i+1,j) && D14(i,j)>D15(i-1,j+1) && D14(i,j)>D15(i,j+1) && D14(i,j)>D15(i+1,j+1))
if(D14(i,j)>threshold)
E13(i,j)=1;
end
end
end
elseif(D14(i,j)<D14(i-1,j-1) && D14(i,j)<D14(i,j-1) && D14(i,j)<D14(i+1,j-1) && D14(i,j)<D14(i-1,j) && D14(i,j)<D14(i+1,j) && D14(i,j)<D14(i-1,j+1) && D14(i,j)<D14(i,j+1) && D14(i,j)<D14(i+1,j+1))
if(D14(i,j)<D13(i,j) && D14(i,j)<D13(i-1,j-1) && D14(i,j)<D13(i,j-1) && D14(i,j)<D13(i+1,j-1) && D14(i,j)<D13(i-1,j) && D14(i,j)<D13(i+1,j) && D14(i,j)<D13(i-1,j+1) && D14(i,j)<D13(i,j+1) && D14(i,j)<D13(i+1,j+1))
if(D14(i,j)<D15(i,j) && D14(i,j)<D15(i-1,j-1) && D14(i,j)<D15(i,j-1) && D14(i,j)<D15(i+1,j-1) && D14(i,j)<D15(i-1,j) && D14(i,j)<D15(i+1,j) && D14(i,j)<D15(i-1,j+1) && D14(i,j)<D15(i,j+1) && D14(i,j)<D15(i+1,j+1))
if(D14(i,j)>threshold)
E13(i,j)=-1;
end
end
end
end
end
end
这样就得到了初步得到了第一层的候选特征点
figure
>> imshow(D11,[0 1])
>> figure
>> imshow(D12,[0 1])
>> figure
>> imshow(D13,[0 1])
>> figure
>> imshow(D14,[0 1])
>> figure
>> imshow(D15,[0 1])
这是D11
这是D12
这是D13
这是D14
这是D15
%下面显示E11 即极大极小值 特征点
figure
>> imshow(E11,[-1 1])
>> figure
>> imshow(E12,[-1 1])
>> figure
>> imshow(E13,[-1 1])
白色表示极大值 黑色表示极小值
但是E12 E13这两幅 出来就看不到特征点了 我在想 是不是 一层都只有一幅有特征点的图出来 后面的两幅图都看不到特征点
下面是E12 和 E13
可以看到
两幅都没有特征点出来 好奇怪明明都是26个像素比较 容我想想 不过既然只有E11有特征点出来 我们就暂时第一层就用E11吧 至少出来了一个
第二层 第三层跟上面步骤道理一样：
m2=round(m1/2);
n2=round(n1/2);
C20=imresize(C11,[m2 n2]);
h21=fspecial('gaussian',[5 5],k*0.3);
C21=imfilter(C20,h21,'conv');
h22=fspecial('gaussian',[5 5],k*0.4);
C22=imfilter(C20,h22,'conv');
h23=fspecial('gaussian',[5 5],k*0.5);
C23=imfilter(C20,h23,'conv');
h24=fspecial('gaussian',[5 5],k*0.6);
C24=imfilter(C20,h24,'conv');
h25=fspecial('gaussian',[5 5],k*0.7);
C25=imfilter(C20,h25,'conv');
h26=fspecial('gaussian',[5 5],k*0.8);
C26=imfilter(C20,h26,'conv');
m3=round(m2/2);
n3=round(n2/2);
C30=imresize(C21,[m3 n3]);
h31=fspecial('gaussian',[5 5],k^2*0.3);
C31=imfilter(C30,h31,'conv');
h32=fspecial('gaussian',[5 5],k^2*0.4);
C32=imfilter(C30,h32,'conv');
h33=fspecial('gaussian',[5 5],k^2*0.5);
C33=imfilter(C30,h33,'conv');
h34=fspecial('gaussian',[5 5],k^2*0.6);
C34=imfilter(C30,h34,'conv');
h35=fspecial('gaussian',[5 5],k^2*0.7);
C35=imfilter(C30,h35,'conv');
h36=fspecial('gaussian',[5 5],k^2*0.8);
C36=imfilter(C30,h36,'conv');
>> D21=C21-C22;
D22=C22-C23;
D23=C23-C24;
D24=C24-C25;
D25=C25-C26;
D31=C31-C32;
D32=C32-C33;
D33=C33-C34;
D34=C34-C35;
D35=C35-C36;
>> figure
>> imshow(D21,[0 1])
这是第二层的第一幅DOG高斯尺度空间的图 依次可以显示第二层第三层的所有图。。。
比较第二层的极值点
E21=zeros(m2,n2);
for i=2:m2-1
for j=2:n2-1
if(D22(i,j)>D22(i-1,j-1) && D22(i,j)>D22(i,j-1) && D22(i,j)>D22(i+1,j-1) && D22(i,j)>D22(i-1,j) && D22(i,j)>D22(i+1,j) && D22(i,j)>D22(i-1,j+1) && D22(i,j)>D22(i,j+1) && D22(i,j)>D22(i+1,j+1))
if(D22(i,j)>D21(i,j) && D22(i,j)>D21(i-1,j-1) && D22(i,j)>D21(i,j-1) && D22(i,j)>D21(i+1,j-1) && D22(i,j)>D21(i-1,j) && D22(i,j)>D21(i+1,j) && D22(i,j)>D21(i-1,j+1) && D22(i,j)>D21(i,j+1) && D22(i,j)>D21(i+1,j+1))
if(D22(i,j)>D23(i,j) && D22(i,j)>D23(i-1,j-1) && D22(i,j)>D23(i,j-1) && D22(i,j)>D23(i+1,j-1) && D22(i,j)>D23(i-1,j) && D22(i,j)>D23(i+1,j) && D22(i,j)>D23(i-1,j+1) && D22(i,j)>D23(i,j+1) && D22(i,j)>D23(i+1,j+1))
if(D22(i,j)>threshold)
E21(i,j)=1;
end
end
end
elseif(D22(i,j)<D22(i-1,j-1) && D22(i,j)<D22(i,j-1) && D22(i,j)<D22(i+1,j-1) && D22(i,j)<D22(i-1,j) && D22(i,j)<D22(i+1,j) && D22(i,j)<D22(i-1,j+1) && D22(i,j)<D22(i,j+1) && D22(i,j)<D22(i+1,j+1))
if(D22(i,j)<D21(i,j) && D22(i,j)<D21(i-1,j-1) && D22(i,j)<D21(i,j-1) && D22(i,j)<D21(i+1,j-1) && D22(i,j)<D21(i-1,j) && D22(i,j)<D21(i+1,j) && D22(i,j)<D21(i-1,j+1) && D22(i,j)<D21(i,j+1) && D22(i,j)<D21(i+1,j+1))
if(D22(i,j)<D23(i,j) && D22(i,j)<D23(i-1,j-1) && D22(i,j)<D23(i,j-1) && D22(i,j)<D23(i+1,j-1) && D22(i,j)<D23(i-1,j) && D22(i,j)<D23(i+1,j) && D22(i,j)<D23(i-1,j+1) && D22(i,j)<D23(i,j+1) && D22(i,j)<D23(i+1,j+1))
if(D22(i,j)<-threshold)
E21(i,j)=-1;
end
end
end
end
end
end
E22=zeros(m2,n2);
for i=2:m2-1
for j=2:n2-1
if(D23(i,j)>D23(i-1,j-1) && D23(i,j)>D23(i,j-1) && D23(i,j)>D23(i+1,j-1) && D23(i,j)>D23(i-1,j) && D23(i,j)>D23(i+1,j) && D23(i,j)>D23(i-1,j+1) && D23(i,j)>D23(i,j+1) && D23(i,j)>D23(i+1,j+1))
if(D23(i,j)>D22(i,j) && D23(i,j)>D22(i-1,j-1) && D23(i,j)>D22(i,j-1) && D23(i,j)>D22(i+1,j-1) && D23(i,j)>D22(i-1,j) && D23(i,j)>D22(i+1,j) && D23(i,j)>D22(i-1,j+1) && D23(i,j)>D22(i,j+1) && D23(i,j)>D22(i+1,j+1))
if(D23(i,j)>D24(i,j) && D23(i,j)>D24(i-1,j-1) && D23(i,j)>D24(i,j-1) && D23(i,j)>D24(i+1,j-1) && D23(i,j)>D24(i-1,j) && D23(i,j)>D24(i+1,j) && D23(i,j)>D24(i-1,j+1) && D23(i,j)>D24(i,j+1) && D23(i,j)>D24(i+1,j+1))
if(D23(i,j)>threshold)
E22(i,j)=1;
end
end
end
elseif(D23(i,j)<D23(i-1,j-1) && D23(i,j)<D23(i,j-1) && D23(i,j)<D23(i+1,j-1) && D23(i,j)<D23(i-1,j) && D23(i,j)<D23(i+1,j) && D23(i,j)<D23(i-1,j+1) && D23(i,j)<D23(i,j+1) && D23(i,j)<D23(i+1,j+1))
if(D23(i,j)<D22(i,j) && D23(i,j)<D22(i-1,j-1) && D23(i,j)<D22(i,j-1) && D23(i,j)<D22(i+1,j-1) && D23(i,j)<D22(i-1,j) && D23(i,j)<D22(i+1,j) && D23(i,j)<D22(i-1,j+1) && D23(i,j)<D22(i,j+1) && D23(i,j)<D22(i+1,j+1))
if(D23(i,j)<D24(i,j) && D23(i,j)<D24(i-1,j-1) && D23(i,j)<D24(i,j-1) && D23(i,j)<D24(i+1,j-1) && D23(i,j)<D24(i-1,j) && D23(i,j)<D24(i+1,j) && D23(i,j)<D24(i-1,j+1) && D23(i,j)<D24(i,j+1) && D23(i,j)<D24(i+1,j+1))
if(D23(i,j)<-threshold)
E22(i,j)=-1;
end
end
end
end
end
end
E23=zeros(m2,n2);
for i=2:m2-1
for j=2:n2-1
if(D24(i,j)>D24(i-1,j-1) && D24(i,j)>D24(i,j-1) && D24(i,j)>D24(i+1,j-1) && D24(i,j)>D24(i-1,j) && D24(i,j)>D24(i+1,j) && D24(i,j)>D24(i-1,j+1) && D24(i,j)>D24(i,j+1) && D24(i,j)>D24(i+1,j+1))
if(D24(i,j)>D23(i,j) && D24(i,j)>D23(i-1,j-1) && D24(i,j)>D23(i,j-1) && D24(i,j)>D23(i+1,j-1) && D24(i,j)>D23(i-1,j) && D24(i,j)>D23(i+1,j) && D24(i,j)>D23(i-1,j+1) && D24(i,j)>D23(i,j+1) && D24(i,j)>D23(i+1,j+1))
if(D24(i,j)>D25(i,j) && D24(i,j)>D25(i-1,j-1) && D24(i,j)>D25(i,j-1) && D24(i,j)>D25(i+1,j-1) && D24(i,j)>D25(i-1,j) && D24(i,j)>D25(i+1,j) && D24(i,j)>D25(i-1,j+1) && D24(i,j)>D25(i,j+1) && D24(i,j)>D25(i+1,j+1))
if(D24(i,j)>threshold)
E23(i,j)=1;
end
end
end
elseif(D24(i,j)<D24(i-1,j-1) && D24(i,j)<D24(i,j-1) && D24(i,j)<D24(i+1,j-1) && D24(i,j)<D24(i-1,j) && D24(i,j)<D24(i+1,j) && D24(i,j)<D24(i-1,j+1) && D24(i,j)<D24(i,j+1) && D24(i,j)<D24(i+1,j+1))
if(D24(i,j)<D23(i,j) && D24(i,j)<D23(i-1,j-1) && D24(i,j)<D23(i,j-1) && D24(i,j)<D23(i+1,j-1) && D24(i,j)<D23(i-1,j) && D24(i,j)<D23(i+1,j) && D24(i,j)<D23(i-1,j+1) && D24(i,j)<D23(i,j+1) && D24(i,j)<D23(i+1,j+1))
if(D24(i,j)<D25(i,j) && D24(i,j)<D25(i-1,j-1) && D24(i,j)<D25(i,j-1) && D24(i,j)<D25(i+1,j-1) && D24(i,j)<D25(i-1,j) && D24(i,j)<D25(i+1,j) && D24(i,j)<D25(i-1,j+1) && D24(i,j)<D25(i,j+1) && D24(i,j)<D25(i+1,j+1))
if(D24(i,j)<-threshold)
E23(i,j)=-1;
end
end
end
end
end
end
下面显示E21
figure
>> imshow(E21,[-1 1])
不知你们看到这幅图的四五个黑白点没 就是极值点 我自己也是醉了
后面第三层的就不显示了 我估计就是一片灰色 连这四五个点都不会给我显示出来了。。。
经过上面的分析以及显示的图像 我们知道其实只需要第一层就可以得到候选特征点 所以我只要E11 完整的程序如下所示：
A=imread('F:\orl_zhifangtu\s1.jpg');
T=40;
[A,R]=kuaisusaomiao(A,T);
[m,n,H]=size(A);
C=[];
x=n/2;
y=m/2;
for u=1:m
for v=1:n
i=u;
j=round(sqrt(R^2-(y-u)^2)*(v-x)/R+x);
if(R^2-(y-u)^2<0)
continue;
end
C(u,v,1)=A(i,j,1);
C(u,v,2)=A(i,j,2);
C(u,v,3)=A(i,j,3);
end
end
C=uint8(C);
C=rgb2gray(C);
H=C;
C=double(C);
[m1,n1,h1]=size(C);
k=2^(1/3);
threshold=3;
h11=fspecial('gaussian',[5 5],0.3);
C11=imfilter(C,h11,'conv');
h12=fspecial('gaussian',[5 5],0.4);
C12=imfilter(C,h12,'conv');
h13=fspecial('gaussian',[5 5],0.5);
C13=imfilter(C,h13,'conv');
h14=fspecial('gaussian',[5 5],0.6);
C14=imfilter(C,h14,'conv');
h15=fspecial('gaussian',[5 5],0.7);
C15=imfilter(C,h15,'conv');
h16=fspecial('gaussian',[5 5],0.8);
C16=imfilter(C,h16,'conv');
D11=C11-C12;
D12=C12-C13;
D13=C13-C14;
D14=C14-C15;
D15=C15-C16;
E11=zeros(m1,n1);
for i=3:m1-2
for j=3:n1-2
if(D12(i,j)>D12(i-1,j-1) && D12(i,j)>D12(i,j-1) && D12(i,j)>D12(i+1,j-1) && D12(i,j)>D12(i-1,j) && D12(i,j)>D12(i+1,j) && D12(i,j)>D12(i-1,j+1) && D12(i,j)>D12(i,j+1) && D12(i,j)>D12(i+1,j+1))
if(D12(i,j)>D11(i,j) && D12(i,j)>D11(i-1,j-1) && D12(i,j)>D11(i,j-1) && D12(i,j)>D11(i+1,j-1) && D12(i,j)>D11(i-1,j) && D12(i,j)>D11(i+1,j) && D12(i,j)>D11(i-1,j+1) && D12(i,j)>D11(i,j+1) && D12(i,j)>D11(i+1,j+1))
if(D12(i,j)>D13(i,j) && D12(i,j)>D13(i-1,j-1) && D12(i,j)>D13(i,j-1) && D12(i,j)>D13(i+1,j-1) && D12(i,j)>D13(i-1,j) && D12(i,j)>D13(i+1,j) && D12(i,j)>D13(i-1,j+1) && D12(i,j)>D13(i,j+1) && D12(i,j)>D13(i+1,j+1))
if(D12(i,j)>threshold)
E11(i,j)=1;
end
end
end
elseif(D12(i,j)<D12(i-1,j-1) && D12(i,j)<D12(i,j-1) && D12(i,j)<D12(i+1,j-1) && D12(i,j)<D12(i-1,j) && D12(i,j)<D12(i+1,j) && D12(i,j)<D12(i-1,j+1) && D12(i,j)<D12(i,j+1) && D12(i,j)<D12(i+1,j+1))
if(D12(i,j)<D11(i,j) && D12(i,j)<D11(i-1,j-1) && D12(i,j)<D11(i,j-1) && D12(i,j)<D11(i+1,j-1) && D12(i,j)<D11(i-1,j) && D12(i,j)<D11(i+1,j) && D12(i,j)<D11(i-1,j+1) && D12(i,j)<D11(i,j+1) && D12(i,j)<D11(i+1,j+1))
if(D12(i,j)<D13(i,j) && D12(i,j)<D13(i-1,j-1) && D12(i,j)<D13(i,j-1) && D12(i,j)<D13(i+1,j-1) && D12(i,j)<D13(i-1,j) && D12(i,j)<D13(i+1,j) && D12(i,j)<D13(i-1,j+1) && D12(i,j)<D13(i,j+1) && D12(i,j)<D13(i+1,j+1))
if(D12(i,j)<-threshold)
E11(i,j)=-1;
end
end
end
end
end
end
%下面这部分是把特征点显示在一个原图中 而不是像E11一样显示在一个空矩阵中
[X11,Y11]=find(E11==1 | E11==-1);
imshow(H)
hold on;
plot(Y11,X11,'rx');（这里敲键盘的时候千万别敲错啊 我之前敲错了 换了X11 Y11的位置 结果出来是个亮瞎我铝合金狗眼的图 我自己想了两天 还以为是我哪里编错了）
这就是以上SIFT程序得到的特征点结果 勉强可以吧
后来搜了http://blog.csdn.net/abcjennifer/article/details/7639681她写的MATLAB程序 我真的有一巴掌拍死自己的冲动 人家写得多么简洁 而我的。。。让我自己一刀了断。。。所以你们还是直接下她的程序比较好 短小精悍 又容易懂
我写这个就当是自己的学习笔记好了。。。用来以后 自己边吐血边回忆。。。
内容来自用户分享和网络整理，不保证内容的准确性，如有侵权内容，可联系管理员处理 
相关文章推荐
- Hadoop概述
- 深入Java注解
- 浏览器缓存机制
- LeetCode 180：Consecutive Numbers
- javascript 全选 反选 js代码
- Mysql时间存储类型优缺点？DATETIME?TIMESTAMP?INT?
- 不同机子tcp报文转发
- mac下eclipse重置workspace
- 文科生也能看懂的iptables教程(推荐)
- Python学习igraph
- 回溯
- 【Unity3d】unity接入sharesdk分享本地图片失败记录
- ROS下安装xtion pro live
- winform中使用ReportViewer的时候,找不到报表数据面板.
- 在公网(internet)上建立website时不能用http访问
- 从两张Excel表所想到的
- DES原理
- 数组中出现次数超过一半的数字
- Kalman滤波各参数的一些调整
- 第十四周项目二——用文件保存的学生名单