✅作者简介：热爱科研的Matlab仿真开发者，修心和技术同步精进，matlab项目合作可私信。
🍎个人主页：Matlab科研工作室
🍊个人信条：格物致知。
更多Matlab仿真内容点击👇
⛄ 内容介绍
从场景的主动照明中捕获低光照水平的深度和反射率图像具有广泛的应用。传统上，即使使用对单个光子敏感的检测器，每个像素也需要数百个光子检测来减轻泊松噪声。我们开发了一种稳健的方法来估计深度和反射率，该方法使用每个像素的固定停留时间以及在场景中平均每个像素检测到一个光子的量级。我们的计算图像形成方法将物理上准确的单光子计数统计与现实世界反射率和 3-D 结构中存在的空间相关性的利用相结合。在强背景光下进行的实验表明，我们的方法能够准确地恢复场景深度和反射率，而基于最大似然 (ML) 估计或其近似的传统成像方法会导致噪声较大的图像。对于深度，性能优于信号无关的噪声去除算法，例如应用于像素级 ML 估计的中值滤波或块匹配和 3-D 滤波 (BM3D)；对于反射率，性能类似于信号相关的噪声去除算法，例如泊松非局部稀疏 PCA 和具有方差稳定变换的 BM3D。我们的框架将光子效率提高了传统处理的 100 倍，并且在光栅扫描操作中总采集时间限制下的第一光子成像也有所改进。因此，我们的新成像器将可用于快速、低功耗和抗噪声的主动光学成像，
⛄ 部分代码
%%
% Sample code and data for
% Photon-efficient computational 3D and reflectivity
% imaging with single-photon detectors
% by D.Shin, A.Kirmani, V.Goyal, and J.H.Shapiro, IEEE TCI 2015
% Uses variants of SPIRAL-TAP software developed in the paper
% This is SPIRAL-TAP: Sparse Poisson intensity reconstruction algorithms?theory and practice
% by Z.Harmany, R.Marcia, and R.Willett, IEEE TIP 2012
% When using this software in your publications, please cite both papers..
%%
clc; clear; close all;
load('g1x1_10.mat');
% load('matlab1.mat')
% load('rjmcmc1.mat')
[m,n,T] = size(mm);
Ts = cell(m,n);
for i = 1:m
for j = 1:n
Ts{i,j} = find(mm(i,j,:));
end
end
[nr,nc] = size(Ts);
C0 = zeros(nr,nc);
D0 = zeros(nr,nc);
for i=1:nr
for j=1:nc
ts = Ts{i,j};
C0(i,j) = length(ts);
D0(i,j) = mean(ts);
end
end
unobserved = isnan(D0);
D0(unobserved) = 0;
addpath(genpath([pwd '/fcns']));
AT = @(x) x; A = @(x) x;
C1 = SPIRALTAP_BIN(C0,A,0.5,100, ...
'noisetype','binomial', ... % approx for binomial
'penalty','tv', ...
'maxiter',60,... % `120
'Initialization',C0,...
'AT',AT,...
'monotone',1,...
'miniter',1,...
'stopcriterion',3,...
'tolerance',1e-12,...
'alphainit',0.1,...
'alphaaccept',1e80,...
'logepsilon',1e-10,...
'saveobjective',1,...
'savereconerror',1,...
'savecputime',1,...
'savesolutionpath',0,...
'truth',C0,...
'verbose',0);
max_count = 20;
scales = C1/max_count;
scales(scales>1) = 1;
img_filt = get_rom(Ts);
D_mean = get_thres(Ts,img_filt,scales);
D_mean(isnan(D_mean)) = 0;
AT = @(x) x; A = @(x) x;
D1 = SPIRALTAP_blank(D_mean,A,10, ...
'noisetype','gaussian', ...
'penalty','tv', ...
'maxiter',40,...
'Initialization',img_filt,...
'AT',AT,...
'monotone',1,...
'miniter',1,...
'stopcriterion',3,...
'tolerance',1e-12,...
'alphainit',0.01,...
'alphaaccept',1e80,...
'logepsilon',1e-10,...
'saveobjective',1,...
'savereconerror',1,...
'savecputime',1,...
'savesolutionpath',0,...
'truth',zeros(nr,nc),...
'verbose',0);
range_C = [0,max_count];
range_D = [1350,1550];
% for k = 1:m
% for l = 1:n
% if dep1(k,l)==4096
% dep1(k,l) = 0;
% D0(k,l) = 0;
% D1(k,l) = 0;
% depth_im(k,l) = 0;
% end
% end
% end
figure;
subplot(221); imagesc(C0, range_C);axis image;axis off;colorbar;colormap jet;
title('reflectivity - ML');
subplot(222); imagesc(D0, range_D);axis image;axis off;colorbar;colormap jet;
title('depth - ML');
subplot(223); imagesc(C1, range_C);axis image;axis off;colorbar;colormap jet;
title('reflectivity - proposed');
subplot(224); imagesc(D1, range_D);axis image;axis off;colorbar;colormap jet;
title('depth - proposed');
% rmse1 = rmse(D0,dep1);
% rmse2 = rmse(D1,dep1);
% rmse3 = rmse(depth_im,dep1);
% subplot(222);imagesc(D1, [0,150]);axis image;axis off;colorbar;colormap jet;
% subplot(223);imagesc(C0, [0,10]);axis image;axis off;colorbar;colormap jet;
% subplot(221);imagesc(D0, [0,150]);axis image;axis off;colorbar;colormap jet;
% subplot(224);imagesc(C1, [0,10]);axis image;axis off;colorbar;colormap jet;
%colormap('gray');
⛄ 运行结果
⛄ 参考文献
