一、PSF生成方法
1. 标准PSF生成
% 高斯PSF生成
gauss_psf = fspecial('gaussian', [7 7], 10); % 7x7高斯核,标准差10% 运动模糊PSF生成
motion_psf = fspecial('motion', 21, 45); % 长度21,角度45度% 点扩散函数可视化
figure;
subplot(121), imshow(gauss_psf, []), title('高斯PSF');
subplot(122), imshow(motion_psf, []), title('运动模糊PSF');
2. 自定义PSF生成
% 生成圆形光瞳PSF
N = 256; % 采样点数
[x,y] = meshgrid(-N/2:N/2-1);
r = sqrt(x.^2 + y.^2);
psf = double(r <= N/4); % 半径为N/4的圆形PSF
psf = psf / sum(psf(:)); % 归一化% 显示PSF
figure, imshow(psf, []), title('圆形光瞳PSF');
二、PSF应用示例
1. 图像模糊模拟
% 读取图像
img = imread('cameraman.tif');
img_gray = rgb2gray(img);% 添加高斯模糊
blurred_img = imfilter(img_gray, gauss_psf, 'conv', 'circular');% 添加噪声
noise_level = 0.01;
noisy_img = imnoise(blurred_img, 'gaussian', 0, noise_level^2);% 显示结果
figure;
subplot(131), imshow(img_gray), title('原始图像');
subplot(132), imshow(blurred_img, []), title('模糊图像');
subplot(133), imshow(noisy_img, []), title('含噪模糊图像');
2. 运动模糊模拟
% 参数设置
theta = 30; % 运动角度
length = 21; % 模糊长度% 生成运动PSF
motion_psf = fspecial('motion', length, theta);% 应用运动模糊
motion_blurred = imfilter(img_gray, motion_psf, 'conv', 'circular');% 显示结果
figure, imshowpair(img_gray, motion_blurred, 'montage');
title('运动模糊效果对比');
三、PSF估计与优化
1. 盲去卷积PSF估计
% 初始化参数
init_psf = ones(size(gauss_psf)); % 初始PSF估计
dampar = 0.1*std(noisy_img(:)); % 阻尼参数
weight = ones(size(noisy_img)); % 权重矩阵% 盲去卷积迭代
[restored_img, est_psf] = deconvblind(noisy_img, init_psf, 15, dampar, weight);% 显示结果
figure;
subplot(131), imshow(noisy_img), title('含噪模糊图像');
subplot(132), imshow(est_psf, []), title('估计PSF');
subplot(133), imshow(restored_img), title('恢复图像');
2. 维纳滤波去模糊
% 计算信噪比
snr = 10*log10(var(img_gray(:))/var(noisy_img(:)-img_gray(:)));% 维纳滤波恢复
wiener_restored = deconvwnr(noisy_img, gauss_psf, snr);% 显示结果
figure, imshowpair(noisy_img, wiener_restored, 'montage');
title('维纳滤波恢复效果');
四、高级PSF建模
1. 非各向同性3D PSF
function psf = nonIsotropicPSF(sigma_xyz)% 生成三维非各向同性高斯PSF[X,Y,Z] = ndgrid(linspace(-1,1,64));psf = exp(-(X.^2/(2*sigma_xyz(1)^2) + Y.^2/(2*sigma_xyz(2)^2) + ...Z.^2/(2*sigma_xyz(3)^2)));psf = psf / sum(psf(:));
end% 使用示例
sigma = [0.1, 0.2, 0.15];
psf_3d = nonIsotropicPSF(sigma);
2. 散射介质PSF模拟
function psf = generateScatteringPSF(N, p)% 生成散射介质PSFhadamard = 1;Nin = 2^hadamard;MyMap = 'hot';% 生成随机相位屏phase_screen = exp(1i*2*pi*rand(N,N,Nin));pupil = zeros(N,N);for j=1:Nfor k=1:Nif (j-0.5-N/2)^2 + (k-0.5-N/2)^2 <= (N/2)^2pupil(j,k) = 1;endendend% 计算散斑speckle = zeros(2^p*N, 2^p*N, Nin);for j=1:Ninspeckle(:,:,j) = fft2(padarray(bruit(:,:,j), [2^(p)*(N/2)-N/2, 2^(p)*(N/2)-N/2]));endPSF = intensity_speckle;
end
五、PSF性能验证
1. 艾里斑验证
% 理想艾里斑参数
lambda = 0.6328e-6; % 波长
f = 100e-3; % 焦距
D = 10e-3; % 孔径直径
v = (pi*D/lambda/f) * (1024/2); % 空间频率% 解析解计算
airy_disk = (2*besselj(1,v)./v).^2;% MATLAB数值仿真
pupil = double(imcircle(1024, 512)); % 创建圆形孔径
fft_pupil = fftshift(fft2(pupil));
simulated_airy = abs(fft_pupil).^2;
simulated_airy = simulated_airy / max(simulated_airy(:));% 对比分析
figure;
subplot(121), imshow(airy_disk, []), title('解析解艾里斑');
subplot(122), imshow(simulated_airy, []), title('数值仿真艾里斑');
2. PSNR/SSIM评估
% 恢复质量评估
original = imread('lena.png');
restored = deconvwnr(noisy_img, gauss_psf);psnr_value = psnr(restored, original);
ssim_value = ssim(restored, original);disp(['PSNR: ', num2str(psnr_value), ' dB']);
disp(['SSIM: ', num2str(ssim_value)]);
参考代码 matlab点扩散函数程序 www.youwenfan.com/contentcnq/63743.html
六、工程优化技巧
-
PSF尺寸选择:PSF尺寸应至少为图像尺寸的1/10
-
噪声抑制:使用
dampar参数控制噪声放大(建议值=0.1*噪声标准差) -
边界处理:采用
edgetaper函数减少边缘效应tapered_img = edgetaper(noisy_img, gauss_psf); -
多通道处理:对彩色图像分通道处理
for c = 1:3[restored(:,:,c), psf(:,:,c)] = deconvblind(noisy_img(:,:,c), init_psf); end
七、扩展应用
-
显微成像系统:结合Zernike多项式模拟像差
zernike = zernike(4, 50); % 4阶球差 aberrated_psf = applyZernike(psf, zernike); -
深度学习辅助:训练PSF估计网络
net = alexnet; net.Layers(end-2) = convolution2dLayer([3 3], 1, 'Padding', 'same'); net = trainNetwork(data, labels, net);