基于机器学习的虹膜识别系统设计
设计的虹膜识别系统流程图如图 1 所示，在图像的预处理过程中主要包括虹膜定位、虹膜区域提取、虹膜区域极坐标变换和归一化处理。最后采用SVM识别方法实现虹膜识别。
图1 虹膜识别系统流程图
虹膜识别主要包括虹膜定位、虹膜裁剪、极坐标变换、归一化、HOG特征提取、虹膜识别六个主要步骤。本系统采用MATLAB作为开发工具实现
图 2 虹膜识别系统设计图
1.虹膜定位
1.1霍夫变换——虹膜快速定位
虹膜定位的目的是较为精确地估计两个圆的中心和半径。大多数经典算法是根据瞳孔的灰度值较低，从而容易找到一个阈值来对图像二值化，进而找到虹膜的内边界，但是在遇到较强烈的光照情况下，瞳孔会有斑点，使阈值选择变得不太容易。而霍夫变换是一种基于参数的投票机制算法，它通过判断图像空间中的每一个边界点是否满足参数空间构成的可能轨迹，以一种累加投票的方式求得边界参数，并根据参数方程将满足参数空间的边界点连接起来，达到目标轮廓跟踪拟合的目的。投票机制是针对所有检测到的边缘点，计算量大是霍夫变换的一个严重缺点。霍夫变换在计算机视觉中广泛应用于轮廓跟踪和图形拟合领域。从理论上来说，霍夫变换可以用来检测任意已知形状的目标。本文采用基于霍夫变换实现虹膜定位，并改进该方法，提高虹膜定位的速度。霍夫变换用于虹膜定位的过程为∶首先，对原虹膜图像用一个边缘检测算子检测其边缘，这里采用经典的Canny 算子，并且得到一个二值边界图像I（xi，yi），（xi，yi）为所有边界点的位置，其中i=1，2，n，n为边界点的个数。然后通过一种投票累加机制用来估计圆心（xj，vj）和半径r。
霍夫变换的定义为：
其中，(xj,yj,r)为圆心坐标和半径大小,构成参数组对应累加器，其值的大小表示该参数组得票多少。得票的多少取决于有多少个边界点落在由该参数组构成的圆上。如果某一边界点(xi,yi,)落在了参数组(xj,yj,r)对应的圆上，那么对应的累加器H(xj,yj,r)的票数就会自动加1。式判断边界点(xi, yi)是否落在(xj,yj ,r)圆上：
其中，
最后统计H (xi, yi , xj,yj ,r)的票数，得票最多的H (xi, yi , xj,yj ,r)所对应的(xj,yj,r)即所求圆的圆心和半径。
霍夫变换的主要缺点在于其计算时间比较长。在霍夫变换中，投票机制要考虑图像中所有边缘点(xi,yi,)到可能的圆心(x0,y0)和半径r的值。霍夫变换的另一个缺点是如果有太多的边缘点，使得搜索范围增大。本文提出一种基于感兴趣区域的虹膜定位算法，该算法改变了以往算法先定位虹膜内边界，再定位外边界的传统思路，而是先定位虹膜的外边界。在搜索边界之前，先定义感兴趣区域，使得搜索范围变小，再去掉一些干扰点，进一步缩小了搜索的范围，有效地解决了霍夫变换运行时间长的问题。
1.2外圆定位
具体算法步骤如下。
(1）在虹膜图像中，取一条水平区域HF (Horizon Field)，定义为从虹膜图像的中心行开始向下取四分之一的图像区域,如图所示。
(2）根据直方图方法得到一个阈值T，用于对H进行二值化，如图9-6(b）所示，并且找到二值化后水平区域中的最大连通区域LUF (The Largest Union Filed)，如图所示。
(3）计算原图像中LUF中的像素的平均值a和方差v,考虑原图像中像素值在[a-v,a+v]的所有像素作为虹膜像素，如图所示。
(4）对原图像中，被LUF包含的区域应用 Canny 算子提取图像的边界，如图所示。
(5）经过边缘提取后得到了一些边缘部分，在这些边缘中忽略掉一些小的边缘，如图所示。
(6）最后通过霍夫变换对这些边缘部分进行处理,得到虹膜的外边界,如图所示。
1.3 内圆定位
(1）计算虹膜区域内部所有像素的灰度平均值 avg，并且以 avg 作为区分虹膜和瞳孔区域的阙值。
(2)根据阙值对图像进行二值化,如图所示。
(3）对图像进行滤波，去除一些噪声，如图所示。通过Canny算子找到瞳孔区域的边界点,如图所示。
(5）对处理后的图像用霍夫变换找到虹膜的内边界，如图所示。
数据库：中科院虹膜数据库CASIA-Iris
一共249类
这里只选取两类进行分析
clc;
clear all;
%% 批量处理数据集
pwd='D:\Desktop\基于机器学习的虹膜识别系统设计';
currentPath = pwd; % 获得当前的工作目录
imds = imageDatastore(fullfile(pwd,'虹膜库/'),'IncludeSubfolders',true,'LabelSource','foldernames'); % 载入所有图片集合
numImages = length(imds.Files); %图片总的张数
for i=1:numImages
eye0=imread(imds.Files{i,1});
Result=preprocess(eye0);
path=strcat('处理后数据\',imds.Files{i,1}(32:35),string(i),'.jpg')
imwrite(Result,path);
end
最后预处理结果
HOG+SVM分类结果
eye0=imread('D:\Desktop\基于机器学习的虹膜识别系统设计\虹膜库\001\L\S1001L02.jpg');
eye1=eye0;
%% 提取眼球
[imgN,imgM]=size(eye1);
for i=1:imgN
for j=1:imgM
if eye1(i,j)<50
eye1(i,j)=0;
else
eye1(i,j)=255;
end
end
end
img2=eye1;
figure;
imshow(img2);
title('提取眼球');
%% 对图像进行滤波
img3=medfilt2(img2);
figure;
imshow(img3);
title('中值滤波');
se=strel('square',3);
img31=imerode(img3,se);
figure;
imshow(img31);
title('形态学处理-膨胀腐蚀');
%% 运用 Canny 算子提取边缘
img4=edge(img31, 'canny');
figure
imshow(img4);
title(' Canny 算子提取边缘')
%% 运用霍夫变换对图像进行圆拟合
[C,HM]=Houghcircle(img4,[30,140]);%%画内圆
I=eye0;
x0=C(1);
yO=C(2);
r=C(3);
DrawCircle(I,C);
%% 画外圆
img=I;
C1(3)=C(3)+55;
C1(1)=C(1);
C1(2)=C(2);
DrawCircle(img,C1);
%% 虹膜裁剪
img4=I;
img5=img4;
for i=1:imgN
for j=1:imgM
if((i-C(2))^2+(j-C(1))^2)>C1(3)^2||((i-C(2))^2+(j-C(1))^2)<C(3)^2
img5(i,j)=255;
end
end
end
figure;
imshow(img5);
title('虹膜裁剪')
%% 极坐标变换
img6=img5;
eyeizb=Expan_normalization(img6,C(1),C(2),C(3),C(1),C(2),C1(3));
figure
imshow(eyeizb)
title('极坐标变换');
%% 虹膜图像规范化
[imgN,imgM]=size(eyeizb);
for i=1:imgN/2
for j=1:imgM
eyegfh(i,j)=eyeizb(i,j);
end
end
figure
imshow(uint8(eyegfh));
title('规范化虹膜32*512');
clc;
clear all;
%% 划分数据集
pwd='D:\Desktop\基于机器学习的虹膜识别系统设计';
currentPath = pwd; % 获得当前的工作目录
imds = imageDatastore(fullfile(pwd,'处理后数据/'),'IncludeSubfolders',true,'LabelSource','foldernames'); % 载入所有图片集合
[imdsTrain,imdsTest]=splitEachLabel(imds,0.8);
N1=length(imdsTrain.Files);
N2=length(imdsTest.Files);
% 初始化特征矩阵
imageSize = [32 512];% 对所有图像进行此尺寸的缩放
I = readimage(imdsTrain,1);
scaleImage = imresize(I,imageSize);
m =size(scaleImage,3);%如果是彩色图片转成灰度图片，灰度图片则不处理
if m>1
scaleImage=rgb2gray(scaleImage);
end
% hog的参数是以下三个，
hog_CellSize=8;
hog_BlockSize=3;
NumBins=9;
% 构造和初始化两种特征向量矩阵
numTrainImages = length(imdsTrain.Files); %训练集样本数
numTestImages = length(imdsTest.Files); %测试集样本数
% featuresTrain 是HOG特征
[hogfeaturestrain, visualization] = extractHOGFeatures(scaleImage,'CellSize',[hog_CellSize hog_CellSize],'BlockSize',[hog_BlockSize hog_BlockSize],'NumBins',NumBins);
featuresTrain = zeros(numTrainImages,size(hogfeaturestrain,2),'double');
% 对每张图片提取两种特征向量
for i = 1:numTrainImages
imageTrain = readimage(imdsTrain,i);
imageTrain = imresize(imageTrain,imageSize); % 缩放
m =size(imageTrain,3);
if m>1
imageTrain=rgb2gray(imageTrain);
end
% HOG
featuresTrain(i,:) = extractHOGFeatures(imageTrain,'CellSize',[hog_CellSize hog_CellSize],'BlockSize',[hog_BlockSize hog_BlockSize],'NumBins',NumBins);
end
% 所有训练图像标签
trainLabels = imdsTrain.Labels; % 训练集的样本标签，是mx1的矩阵，是文字的
double_trainLabels=double(imdsTrain.Labels);% 将标签转成数值
%% 提取测试集的特征和特征标签
imageSize = [32 512];% 对所有图像进行此尺寸的缩放 参数1：图片大小
I = readimage(imdsTest,1);
scaleImage = imresize(I,imageSize);
m =size(scaleImage,3);% 如果是彩色图片转成灰度图片，灰度图片则不处理
if m>1
scaleImage=rgb2gray(scaleImage);
end
numTestImages = length(imdsTest.Files); % 测试集样本数
% featuresTest2是HOG特征
[hogfeaturestest, visualization] = extractHOGFeatures(scaleImage,'CellSize',[hog_CellSize hog_CellSize],'BlockSize',[hog_BlockSize hog_BlockSize],'NumBins',NumBins);
featuresTest2 = zeros(numTestImages,size(hogfeaturestest,2),'double');
for i = 1:numTestImages
imageTest = readimage(imdsTest,i);
imageTest = imresize(imageTest,imageSize); %缩放
m =size(imageTest,3);
if m>1
imageTest=rgb2gray(imageTest);
end
% HOG
featuresTest2(i,:) = extractHOGFeatures(imageTest,'CellSize',[hog_CellSize hog_CellSize],'BlockSize',[hog_BlockSize hog_BlockSize],'NumBins',NumBins);
end
featuresTest = featuresTest2;% 把两种特征横向地拼接起来
% 所有测试图像标签
testLabels = imdsTest.Labels; % 训练集的样本标签每一类的数量，是mx1的矩阵
double_testLabels=double(imdsTest.Labels);% 将标签转成数值
%% svm分类
SVMModel=fitcsvm(featuresTrain,double_trainLabels)
%% 测试
s=0;
predict_testLabels=SVMModel.predict(featuresTest);
for i=1:length(predict_testLabels)
if predict_testLabels(i)==double_testLabels(i)
s=s+1;
end
end
acc=s/length(predict_testLabels);
disp('测试集正确率')
acc
confusionchart(double_testLabels, predict_testLabels);