MATLAB 虹膜识别系统：基于霍夫变换的定位与分割

使用MATLAB实现的虹膜识别系统，利用霍夫变换进行虹膜定位和分割

classdef IrisRecognitionSystem
    %IRISRECOGNITIONSYSTEM 虹膜识别系统类
    %   使用霍夫变换进行虹膜定位和分割
    
    properties
        OriginalImage
        GrayImage
        ProcessedImage
        PupilCircle
        IrisCircle
        NormalizedIris
    end
    
    methods
        function obj = IrisRecognitionSystem()
            %IRISRECOGNITIONSYSTEM 构造此类的实例
            obj.PupilCircle = [];
            obj.IrisCircle = [];
        end
        
        function obj = loadImage(obj, imagePath)
            %LOADIMAGE 加载图像
            obj.OriginalImage = imread(imagePath);
            if size(obj.OriginalImage, 3) == 3
                obj.GrayImage = rgb2gray(obj.OriginalImage);
            else
                obj.GrayImage = obj.OriginalImage;
            end
            figure, imshow(obj.OriginalImage), title('原始图像');
        end
        
        function obj = preprocessImage(obj)
            %PREPROCESSIMAGE 图像预处理
            % 应用中值滤波减少噪声
            obj.ProcessedImage = medfilt2(obj.GrayImage, [5 5]);
            
            % 对比度增强
            obj.ProcessedImage = imadjust(obj.ProcessedImage);
            
            figure, imshow(obj.ProcessedImage), title('预处理后的图像');
        end
        
        function obj = detectPupil(obj)
            %DETECTPUPIL 使用霍夫变换检测瞳孔
            
            % 二值化图像
            threshold = graythresh(obj.ProcessedImage);
            binaryImage = imbinarize(obj.ProcessedImage, threshold * 0.8);
            binaryImage = imcomplement(binaryImage);
            
            % 形态学操作
            se = strel('disk', 3);
            binaryImage = imopen(binaryImage, se);
            binaryImage = imclose(binaryImage, se);
            
            % 填充孔洞
            binaryImage = imfill(binaryImage, 'holes');
            
            % 边缘检测
            edges = edge(binaryImage, 'canny');
            
            % 霍夫变换检测圆
            radii = 20:5:80; % 瞳孔半径范围
            h = circle_hough(edges, radii, 'same', 'normalise');
            
            % 找到最可能的圆
            peaks = circle_houghpeaks(h, radii, 'nhoodxy', 15, 'nhoodr', 15, 'npeaks', 1);
            
            if ~isempty(peaks)
                obj.PupilCircle = peaks;
                fprintf('检测到瞳孔: 中心(%d, %d), 半径%d\n', peaks(1), peaks(2), peaks(3));
                
                % 绘制检测结果
                figure, imshow(obj.OriginalImage);
                hold on;
                plot(peaks(1), peaks(2), 'r+', 'MarkerSize', 20, 'LineWidth', 2);
                theta = 0:0.01:2*pi;
                x = peaks(3) * cos(theta) + peaks(1);
                y = peaks(3) * sin(theta) + peaks(2);
                plot(x, y, 'r-', 'LineWidth', 2);
                title('瞳孔检测结果');
                hold off;
            else
                error('未能检测到瞳孔');
            end
        end
        
        function obj = detectIris(obj)
            %DETECTIRIS 使用霍夫变换检测虹膜边界
            
            if isempty(obj.PupilCircle)
                error('请先检测瞳孔');
            end
            
            % 创建感兴趣区域(ROI)掩码
            [height, width] = size(obj.ProcessedImage);
            [X, Y] = meshgrid(1:width, 1:height);
            
            % 创建以瞳孔为中心的圆形掩码
            pupilX = obj.PupilCircle(1);
            pupilY = obj.PupilCircle(2);
            mask = sqrt((X - pupilX).^2 + (Y - pupilY).^2) <= 150;
            
            % 应用掩码
            roiImage = obj.ProcessedImage;
            roiImage(~mask) = 0;
            
            % 边缘检测
            edges = edge(roiImage, 'canny');
            
            % 霍夫变换检测圆 (虹膜半径通常比瞳孔大)
            radii = obj.PupilCircle(3)+30:5:150;
            h = circle_hough(edges, radii, 'same', 'normalise');
            
            % 找到最可能的圆
            peaks = circle_houghpeaks(h, radii, 'nhoodxy', 15, 'nhoodr', 15, 'npeaks', 1);
            
            if ~isempty(peaks)
                obj.IrisCircle = peaks;
                fprintf('检测到虹膜: 中心(%d, %d), 半径%d\n', peaks(1), peaks(2), peaks(3));
                
                % 绘制检测结果
                figure, imshow(obj.OriginalImage);
                hold on;
                
                % 绘制瞳孔
                plot(obj.PupilCircle(1), obj.PupilCircle(2), 'r+', 'MarkerSize', 20, 'LineWidth', 2);
                theta = 0:0.01:2*pi;
                x_pupil = obj.PupilCircle(3) * cos(theta) + obj.PupilCircle(1);
                y_pupil = obj.PupilCircle(3) * sin(theta) + obj.PupilCircle(2);
                plot(x_pupil, y_pupil, 'r-', 'LineWidth', 2);
                
                % 绘制虹膜
                plot(peaks(1), peaks(2), 'g+', 'MarkerSize', 20, 'LineWidth', 2);
                x_iris = peaks(3) * cos(theta) + peaks(1);
                y_iris = peaks(3) * sin(theta) + peaks(2);
                plot(x_iris, y_iris, 'g-', 'LineWidth', 2);
                
                title('虹膜检测结果');
                hold off;
            else
                error('未能检测到虹膜');
            end
        end
        
        function obj = extractIris(obj)
            %EXTRACTIRIS 提取并归一化虹膜区域
            
            if isempty(obj.PupilCircle) || isempty(obj.IrisCircle)
                error('请先检测瞳孔和虹膜');
            end
            
            % 创建极坐标网格
            [height, width] = size(obj.ProcessedImage);
            [X, Y] = meshgrid(1:width, 1:height);
            
            % 计算每个点到瞳孔中心的距离和角度
            dx = X - obj.PupilCircle(1);
            dy = Y - obj.PupilCircle(2);
            angles = atan2(dy, dx);
            radii = sqrt(dx.^2 + dy.^2);
            
            % 归一化参数
            normalizedWidth = 256;
            normalizedHeight = 64;
            
            % 创建归一化虹膜图像
            obj.NormalizedIris = zeros(normalizedHeight, normalizedWidth);
            
            for y = 1:normalizedHeight
                for x = 1:normalizedWidth
                    % 计算对应的极坐标
                    r = obj.PupilCircle(3) + (y/normalizedHeight) * ...
                        (obj.IrisCircle(3) - obj.PupilCircle(3));
                    theta = (x/normalizedWidth) * 2 * pi;
                    
                    % 转换为笛卡尔坐标
                    x_cart = round(obj.PupilCircle(1) + r * cos(theta));
                    y_cart = round(obj.PupilCircle(2) + r * sin(theta));
                    
                    % 确保坐标在图像范围内
                    if x_cart >= 1 && x_cart <= width && y_cart >= 1 && y_cart <= height
                        obj.NormalizedIris(y, x) = obj.ProcessedImage(y_cart, x_cart);
                    end
                end
            end
            
            % 显示归一化虹膜
            figure, imshow(obj.NormalizedIris), title('归一化虹膜');
            
            % 增强对比度
            obj.NormalizedIris = imadjust(obj.NormalizedIris);
            figure, imshow(obj.NormalizedIris), title('增强后的归一化虹膜');
        end
        
        function features = extractFeatures(obj)
            %EXTRACTFEATURES 从归一化虹膜中提取特征
            
            if isempty(obj.NormalizedIris)
                error('请先提取虹膜');
            end
            
            % 使用Gabor滤波器提取纹理特征
            wavelength = 8;
            orientation = 0:30:150;
            g = gabor(wavelength, orientation);
            
            % 应用Gabor滤波器
            gabormag = imgaborfilt(obj.NormalizedIris, g);
            
            % 显示Gabor滤波结果
            figure;
            for i = 1:length(g)
                subplot(2,3,i), imshow(gabormag(:,:,i), []);
                title(sprintf('方向 %d度', orientation(i)));
            end
            
            % 提取特征向量
            features = [];
            for i = 1:length(g)
                mag = gabormag(:,:,i);
                features = [features, mean(mag(:)), std(mag(:))];
            end
            
            fprintf('提取的特征向量: \n');
            disp(features);
        end
    end
end

辅助函数

为了完整实现虹膜识别系统，您还需要以下辅助函数。这些函数应该保存在MATLAB路径中。

circle_hough.m

function [h, varargout] = circle_hough(b, r, varargin)
%CIRCLE_HOUGH Hough transform for circles
%   [H, R] = CIRCLE_HOUGH(B, R) computes the Hough transform for circles
%   of radii in the vector R in the binary image B.
%
%   [H, R] = CIRCLE_HOUGH(B, R, OPTION1, VALUE1, ...) accepts additional
%   options:
%     'same'      - keep H the same size as B (default)
%     'normalise' - normalise the transform by the radius

% 参数处理
opts = struct('same', false, 'normalise', false);
if nargin > 2
    opts = parse_args(opts, varargin);
end

% 初始化累加器
[h, w] = size(b);
n_r = length(r);
h = zeros(h, w, n_r);

% 找到所有边缘点
[y, x] = find(b);

% 对每个半径进行霍夫变换
for i = 1:n_r
    radius = r(i);
    
    % 创建圆模板
    angles = 0:0.1:2*pi;
    template_x = round(radius * cos(angles));
    template_y = round(radius * sin(angles));
    
    % 移除重复点
    template = unique([template_x', template_y'], 'rows');
    
    % 对每个边缘点投票
    for j = 1:length(x)
        centers_x = x(j) - template(:,1);
        centers_y = y(j) - template(:,2);
        
        % 移除超出图像边界的点
        valid = centers_x >= 1 & centers_x <= w & centers_y >= 1 & centers_y <= h;
        centers_x = centers_x(valid);
        centers_y = centers_y(valid);
        
        % 投票
        for k = 1:length(centers_x)
            h(centers_y(k), centers_x(k), i) = h(centers_y(k), centers_x(k), i) + 1;
        end
    end
    
    % 归一化
    if opts.normalise
        h(:,:,i) = h(:,:,i) / radius;
    end
end

% 输出半径向量
if nargout > 1
    varargout{1} = r;
end
end

function opts = parse_args(opts, args)
%PARSE_ARGS 解析参数
n = length(args);
for i = 1:2:n
    if isfield(opts, lower(args{i}))
        opts.(lower(args{i})) = args{i+1};
    end
end
end

circle_houghpeaks.m

function peaks = circle_houghpeaks(h, radii, varargin)
%CIRCLE_HOUGHPEAKS Find peaks in circle Hough transform
%   PEAKS = CIRCLE_HOUGHPEAKS(H, RADII, PARAM1, VAL1, ...) finds peaks in the
%   circle Hough transform accumulator H for circles with radii RADII.

% 默认参数
params = struct('nhoodxy', 15, 'nhoodr', 15, 'npeaks', 1, 'threshold', 0.5);

% 解析输入参数
if nargin > 2
    params = parse_params(params, varargin);
end

% 找到最大值的初始候选
max_val = max(h(:));
threshold = params.threshold * max_val;

peaks = [];
for i = 1:params.npeaks
    % 找到当前最大值
    [r, c, s] = ind2sub(size(h), find(h == max_val, 1));
    
    if isempty(r) || h(r, c, s) < threshold
        break;
    end
    
    % 保存峰值
    peaks = [peaks; c, r, radii(s)];
    
    % 抑制邻域
    r_start = max(1, r - params.nhoodxy);
    r_end = min(size(h,1), r + params.nhoodxy);
    c_start = max(1, c - params.nhoodxy);
    c_end = min(size(h,2), c + params.nhoodxy);
    s_start = max(1, s - params.nhoodr);
    s_end = min(size(h,3), s + params.nhoodr);
    
    h(r_start:r_end, c_start:c_end, s_start:s_end) = 0;
    
    % 更新最大值
    max_val = max(h(:));
end
end

function params = parse_params(params, args)
%PARSE_PARAMS 解析参数
n = length(args);
for i = 1:2:n
    if isfield(params, lower(args{i}))
        params.(lower(args{i})) = args{i+1};
    end
end
end

使用示例

% 创建虹膜识别系统实例
irisSystem = IrisRecognitionSystem();

% 加载图像
irisSystem = irisSystem.loadImage('iris_image.jpg');

% 预处理图像
irisSystem = irisSystem.preprocessImage();

% 检测瞳孔
irisSystem = irisSystem.detectPupil();

% 检测虹膜
irisSystem = irisSystem.detectIris();

% 提取并归一化虹膜
irisSystem = irisSystem.extractIris();

% 提取特征
features = irisSystem.extractFeatures();

参考代码 虹膜识别源代码 www.3dddown.com/cna/55118.html

说明

1. 图像预处理：应用中值滤波和对比度增强来提高图像质量。
2. 瞳孔检测：使用霍夫变换检测瞳孔，这是虹膜的内边界。
3. 虹膜检测：在瞳孔周围区域使用霍夫变换检测虹膜的外边界。
4. 虹膜归一化：将环形虹膜区域转换为矩形区域，以便于后续特征提取。
5. 特征提取：使用Gabor滤波器提取虹膜纹理特征。

注意

1. 此代码需要Image Processing Toolbox。
2. 霍夫变换参数可能需要根据具体图像进行调整。
3. 对于不同的虹膜数据库，可能需要调整预处理步骤和参数。
4. 完整的虹膜识别系统还需要添加匹配和识别模块。