计算机视觉:朗伯光度立体法(Lambertian Photometric Stereo)

计算机视觉:朗伯光度立体法(Lambertian Photometric Stereo)

  • 光度立体法简介
  • 朗伯光度立体法算法原理
  • 朗伯光度立体法matlab程序
  • 示例
    • Albedo图
    • Normal图
    • Re_rendered图
  • 参考文献

光度立体法简介

光度立体法,即Photometric Stereo, 最早是由当时在MIT的人工智能实验室的Robert J. Woodham教授在1978年左右提出。他在1979年的论文《Photometric stereo: A reflectance map technique for determining surface orientation from image intensity》,以及1980年的论文《Photometric Method for Determining Surface Orientation from Multiple Images》中比较系统的阐述了整套理论框架。光度立体法可以根据二维纹理信息提取出三维模型,其典型应用是检测物体表面微小变化。

朗伯光度立体法基于Woodham算法,Woodham在论文中提出三个基本假设。在这三个基本假设下完成整套理论框架的推演。这三个基本假设分别是:

  1. 假定相机是无畸变成像,也就是说必须使用远心镜头或者长焦镜头。
  2. 假定每一个光源发射的光束都是平行且均匀的,也就是说必须使用具有均匀强度的远心照明光源,或者使用远距离的点光源代替。
  3. 物体必须具有朗伯(lambertian)反射特性,即它必须以漫反射的方式反射入射光。

朗伯光度立体法的大致思路是:
当相机和目标物体相对位置固定不变时,使用不同方向的光源照射同一目标物体,相机可以拍摄到目标物体带有不同明暗分布的图像(至少需要三张图),再通过求解基于朗伯反射原理的反射方程组,求解目标表面的法向分布或者albedo图。

带有远心镜头的相机必须与被测物体表面垂直安装,在采集多幅图像时,一定要保证相机和物体不被移动。相反,对于采集至少三张的灰度图像,其每次取像的照明方向必须改变(相对于相机)。对于采集的多张图像中的每一幅图,照明方向必须指定Slants和Tilts两个参数角度,其描述了相对于当前场景的光照角度。为了更好的理解这两个参数含义,我们假定光源射出的光束是平行光,镜头是远心镜头,相机垂直于物体表面。
在这里插入图片描述
在这里插入图片描述

朗伯光度立体法算法原理

根据Lambert模型:
I = ρ L ⋅ N \textbf{I} = \rho\ \textbf{L}\cdot\textbf{N} I=ρ LN
式中, ρ \rho ρ 为表面反射率(albedo),其值介于0 - 1之间,反映了物体表面特性; N \textbf{N} N 为表面法线(normal); L \textbf{L} L 为照明方向。将 ρ \rho ρ N \textbf{N} N 合并起来用 G \textbf{G} G来表示,则有:
I = L T G . \textbf{I} = \textbf{L}^T\textbf{G}. I=LTG.
每个像素位置对应的 G \textbf{G} G是三维向量(方向为normal,模长为albedo),因此单个光源无法求解该方程,至少需要三个不同位置的光源,求解得到三个未知量。
由最小二乘法:
min ⁡ G ∣ ∣ I − L T G ∣ ∣ 2 \mathop{\min}\limits_{\textbf{G}} ||\textbf{I}-\textbf{L}^T\textbf{G}||^2 Gmin∣∣ILTG2
可以求得:
G = ( L L T ) − 1 LI \textbf{G} = (\textbf{L}\textbf{L}^T)^{-1}\textbf{L}\textbf{I} G=(LLT)1LI
我们求得 G \textbf{G} G 的模长就是albedo:
ρ = ∣ ∣ G ∣ ∣ \rho = ||\textbf{G}|| ρ=∣∣G∣∣
G \textbf{G} G 归一化后的单位向量矩阵就是normal:
N = G ∣ ∣ G ∣ ∣ \textbf{N} = \frac{\textbf{G}}{||\textbf{G}||} N=∣∣G∣∣G

朗伯光度立体法matlab程序

注:此代码只涉及核心算法,不包含数据的读入,与结果的输出。

function [Albedo, Normal, Re_rendered] = Photometric_Stereo(data)assert(size(data.imgs, 1)==size(data.s, 1), 'Size mismatched!');
num = size(data.imgs, 1);% Get image dimensions
im = data.imgs{1};
[im_h, im_w, ~] = size(im);
% Initialize T, a im_h-by-im_w-by-num matrix, whose (h, w, :) holds 
% the intensities at (h, w) for all p different lightings
im_T = zeros(im_h, im_w, num);
% Initialize im_R, a im_h-by-im_w-by-num matrix
im_R = zeros(im_h, im_w, num);
% Initialize im_G, a im_h-by-im_w-by-num matrix
im_G = zeros(im_h, im_w, num);
% Initialize im_B, a im_h-by-im_w-by-num matrix
im_B = zeros(im_h, im_w, num);% For each image
for idx = 1:numim = data.imgs{idx};imGray = rgb2gray(im);% Loop thru each pixelfor h = 1:im_hfor w = 1:im_w% If in the maskif data.mask(h, w)im_R(h, w, idx) = im(h, w, 1);im_G(h, w, idx) = im(h, w, 2);im_B(h, w, idx) = im(h, w, 3);im_T(h, w, idx) = imGray(h,w);endendend
end% Initialize Normal, a im_h-by-im_w-by-3 matrix
Normal = zeros(im_h, im_w, 3);
% Initialize Albedo, a im_h-by-im_w-by-1 matrix
Albedo = zeros(im_h, im_w, 3);
% Initialize Re_rendered, a im_h-by-im_w-by-3 matrix
Re_rendered = zeros(im_h, im_w, 3);% Loop thru each location
for h = 1:im_hfor w = 1:im_w% If in the maskif data.mask(h, w)%% Normal % LightingsL = data.s;% IntensitiesI = reshape(im_T(h, w, :), [num, 1]);% Dealing with shadows and highlightsfor k = 1:10max_index = find(I==max(I));I(max_index) = [];L(max_index,:) = [];min_index = find(I==min(I));I(min_index) = [];L(min_index,:) = [];end% Solve surface normals and albedoG = (L.'*L)\(L.'*I);if norm(G) ~= 0% Normalize nn = G./norm(G);elsen = [0; 0; 0];end% SaveNormal(h, w, :) = n;%% Albedo% a_RL = data.s;I_R = reshape(im_R(h, w, :), [num, 1]);for k = 1:10max_index = find(I_R==max(I_R));I_R(max_index) = [];L(max_index,:) = [];min_index = find(I_R==min(I_R));I_R(min_index) = [];L(min_index,:) = [];endG_R = (L.'*L)\(L.'*I_R);a_R = norm(G_R);Albedo(h, w, 1) = a_R;% a_GL = data.s;I_G = reshape(im_G(h, w, :), [num, 1]);for k = 1:10max_index = find(I_G==max(I_G));I_G(max_index) = [];L(max_index,:) = [];min_index = find(I_G==min(I_G));I_G(min_index) = [];L(min_index,:) = [];endG_G = (L.'*L)\(L.'*I_G);a_G = norm(G_G);Albedo(h, w, 2) = a_G;% a_BL = data.s;I_B = reshape(im_B(h, w, :), [num, 1]);for k = 1:10max_index = find(I_B==max(I_B));I_B(max_index) = [];L(max_index,:) = [];min_index = find(I_B==min(I_B));I_B(min_index) = [];L(min_index,:) = [];endG_B = (L.'*L)\(L.'*I_B);a_B = norm(G_B);Albedo(h, w, 3) = a_B;%% Re_renderedRe_rendered(h, w, :) = [a_R;a_G;a_B]*dot(n,[0;0;1]);endend
end

完整程序请移步至此链接下载

示例

Albedo图

反照率图:
反照率图

Normal图

以RGB线性编码的法线贴图:
以RGB线性编码的法线贴图

Re_rendered图

在照明方向和观察方向一致时,利用normal和albedo图重新渲染的图片:

在照明方向和观察方向一致时,利用normal和albedo图重新渲染的图片

参考文献

光度立体简介
Halcon 光度立体法(photometric_stereo)详解z

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mzph.cn/news/587755.shtml

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈email:809451989@qq.com,一经查实,立即删除!

相关文章

啊哈c语言——4.10(练习)

1&#xff0e;请尝试用for循环打印下面的图形。 #include <stdio.h> #include <stdlib.h> int main() {int a,b,c,d,e;for(a 1;a < 10;a){if(a < 5){b a * 2 - 1;c 5 - a;}else{b 9 - (a - 5) * 2;c a - 5;}for(d 0;d < c;d ){printf(" "…

【LeetCode每日一题】1154. 一年中的第几天(直接计算+调用库函数)

2023-12-31 文章目录 [1154. 一年中的第几天](https://leetcode.cn/problems/day-of-the-year/)方法一&#xff1a;直接计算思路&#xff1a; 方法二&#xff1a;调用库函数思路 1154. 一年中的第几天 方法一&#xff1a;直接计算 思路&#xff1a; 1.根据所给的字符串&#…

【机器学习】快速入门!关于 Pandas 库的简介和常用方法整理

Pandas Pandas 简介1. 数据加载和存储加载数据&#xff1a;存储数据&#xff1a; 2. 数据清洗3. 数据统计和汇总4. 数据选择和过滤5. 数据合并和连接6. 时间序列处理创建时间序列数据&#xff1a;索引和选择&#xff1a;时间序列分析&#xff1a;时间序列可视化&#xff1a; 7.…

基于SpringBoot实现的前后端分离电影评分项目,功能:注册登录、浏览影片、热门影片、搜索、评分、片单、聊天、动态

一、项目介绍 本项目主要基于SpringBoot、Mybatis-plus、MySQL、Redis实现的影片评分项目。 本系统是前后端分离的&#xff0c;分别由三个子项目构成&#xff1a;java服务端、用户前端、管理员管理前端 关键词&#xff1a;springboot java vue mysql reids websocket 毕业设计…

vue el-select 设置默认值后选项无法切换

vue el-select 设置默认值后选项无法切换 解决方式 change"$forceUpdate()" 添加这个即可 完整的代码 <template><el-form-item label"数据类型"><el-select v-model"queryParams.searchDataType" placeholder"请选择数据…

在Neo4j中实现推荐算法

在Neo4j中实现推荐算法 推荐系统是当今信息过载时代的关键技术&#xff0c;它帮助用户在海量数据中发现对他们可能有用或感兴趣的内容。在社交网络、电子商务和内容平台等多个领域&#xff0c;推荐算法的应用已经变得非常广泛。图数据库如Neo4j因其天然对关系数据的支持&#…

V子型输出一串字符。。。

#include<stdio.h>int chars[100][2] {};//每行要出现的字符编码 int main() {int line;char start;char c, c1;scanf("%d %c",&line,&start);//输出多少行for (int i 0; i < line; i) {c A (start-A i)%26;c1 A (start - A 2*line-i-2) % 26…

【REST2SQL】01RDB关系型数据库REST初设计

0 概念 REST2SQL实现连接数据库&#xff0c;数据库的表或视图即可提供REST的GET\POST\PUT\DELETE请求&#xff0c;SQL可执行SQLECT\INSERT\UPDATE\DELETE语句。 0.1 RDB Relational Database 即关系型数据库&#xff08;简称 RDB&#xff09;是一种以关系&#xff08;即表格…

数据结构,题目笔记

哈希表 线性探测再散列 【算法数据结构&#xff5c;哈希查找&#xff5c;哈希冲突&#xff5c;除留余数法&#xff5c;线形探测法&#xff5c;例题讲解】https://www.bilibili.com/video/BV1514y1P7BK?vd_source1a684a3a1b9d05485b3d6277aeeb705d 【二次探测再散列法】 【【…

Linux sed 命令

你可以使用 sed 命令来在 shell 中找到文件中的特定字符并进行替换。以下是一个简单的示例&#xff1a; 假设要将文件中的所有 old_string 替换为 new_string&#xff0c;可以使用以下命令&#xff1a; sed -i s/dynamic_4d_app\b/dynamic_4d_app_associaqted/g CMakeLists.t…

ESP32入门七(中断)

中断用于处理在程序正常执行期间通过外部事件或者响应软件指令触发时发生的事件。比如&#xff0c;在一段呼吸灯的代码中&#xff0c;正常运行时的结果为LED从暗到亮&#xff0c;再从亮到暗持续地运行。我们可以通过一个中断来控制呼吸灯的运行和停止。使用中断功能&#xff0c…

C Primer Plus (中文版)第11章编程练习 参考答案(仅供参考~)

C Primer Plus &#xff08;中文版&#xff09;第11章编程练习 参考答案(仅供参考~) &#x1f334; C Primer Plus第11章编程练习~ 加油加油&#xff01;&#x1f36d; &#x1f36d;感觉这一章&#xff0c;比较吃力~ 很迟没有更新了&#xff0c;也有自己的原因 ~ 2023的最后一…

OpenCV-11颜色通道的分离与合并

本次我们使用两个比较重要的API split&#xff08;mat&#xff09;将图像的通道进行分割。 merge&#xff08;(ch1&#xff0c;ch2&#xff0c;ch3)&#xff09;将多个通道进行融合。 示例代码如下&#xff1a; import cv2 import numpy as npimg np.zeros((480, 640, 3),…

uni-app App.vue生命周期全局样式全局存储globalData

锋哥原创的uni-app视频教程&#xff1a; 2023版uniapp从入门到上天视频教程(Java后端无废话版)&#xff0c;火爆更新中..._哔哩哔哩_bilibili2023版uniapp从入门到上天视频教程(Java后端无废话版)&#xff0c;火爆更新中...共计23条视频&#xff0c;包括&#xff1a;第1讲 uni…

Docker之网络配置

目录 1.网络概念 网络相关的有ip,子网掩码,网关,DNS,端口号 1.1 ip是什么? ip是唯一定位一台网上计算机 Ip地址的分类: IPV4: 4字节32位整数&#xff0c;并分成4段8位的二进制数&#xff0c;每8位之间用圆点隔开&#xff0c;每8位整数可以转换为一个0~255的十进制整数 【例…

windows server 2022 启用SYN攻击保护

2023.12.28 SYN攻击是什么&#xff1a; SYN攻击是黑客攻击的常用手段&#xff0c;也是最容易被利用的一种攻击手法&#xff0c;属于DDoS攻击的一种。它利用TCP协议缺陷&#xff0c;通过发送大量的半连接请求&#xff0c;耗费CPU和内存资源。 SYN攻击包括大量TCP连接的第一个包&…

Leetcode 第 378 场双周赛 Problem D 回文串重新排列查询(Java + 区间相交关系 + 前缀和)

文章目录 题目思路Java 区间相交关系 前缀和:第 1 步&#xff1a;第 2 步&#xff1a;第 3 步&#xff1a;第 4 步&#xff1a; 复杂度Code 题目 100129. 回文串重新排列查询给你一个长度为 偶数 n &#xff0c;下标从 0 开始的字符串 s 。同时给你一个下标从 0 开始的二维整…

算法题Python常用内置函数、方法、技巧汇总(其五:堆/优先队列)

文章目录 优先队列相关操作堆化入堆出堆获取堆顶元素小根堆与大根堆 华为OD算法/大厂面试高频题算法练习冲刺训练 优先队列相关操作 注意&#xff0c;优先队列&#xff08;priority queue&#xff09;也叫做堆&#xff08;heap&#xff09;。谈到优先队列时&#xff0c;一般强调…

Nginx屏蔽垃圾邮件骚扰IP的方法

原文地址&#xff1a;Nginx屏蔽垃圾邮件骚扰IP的方法 本文介绍了如何下载并引入deny-ips.conf配置文件&#xff0c;以及如何定制403页面&#xff0c;避免误杀合法访问者。 最近&#xff0c;很多人都遭受到垃圾邮件的骚扰&#xff0c;让我们无法正常地观看和回复重要信息。在这种…

Python面向对象-三大特性

一 三大特性-封装 面向对象思想有三大特性&#xff1a;封装、继承、多态。 封装&#xff1a;将属性和方法放到一起做为一个整体&#xff0c;然后通过实例化对象来处理&#xff0c;这样隐藏内部实现细节&#xff0c;只需要和对象及其属性和方法交互就可以了。 为了更好的封装…