CN102779354B - 基于光度立体技术的中医望诊信息表面三维重建方法 - Google Patents

Abstract

本发明属于中医望诊的图像处理领域，特别涉及基于光度立体技术的中医望诊信息三维表面重建领域。本发明基于不同光源方向下的两组图像，即以黑球为拍摄对象的高光黑球图像和以表面待重建对象为拍摄对象的表面待重建对象图像。首先对这两组图像进行预处理，包括图片校正、图片裁剪、以及制作Mask图片，并计算光源方向向量L；然后计算待重建对象三维重建图中任意像素点q的单位法向量；最后计算待重建对象三维重建图的深度信息矩阵，并建立深度信息Z与待重建对象三维重建图中各像素点的对应关系，绘制待重建对象的三维重建图。本方法可实现望诊信息表面的三维重建，进行望诊信息的采集和表达，并利用重建的三维视图多角度观察其表面形态。

Description

基于光度立体技术的中医望诊信息表面三维重建方法

技术领域

本发明属于中医望诊的图像处理领域，特别涉及基于光度立体技术的中医望诊信息三维表面重建领域。

背景技术

舌诊和面诊是我国传统医学四诊“望、闻、问、切”中望诊的重要内容，是中医诊法最重要的特色之一。特别是舌象的表面形态以及舌体上的裂纹和舌乳头的分布大小等细节信息与中医舌诊的某些病变有重要联系，对其诊断和治疗有重要的现实意义。鉴于传统的望诊方法主要依靠医生目测进行判断分析，缺乏客观评价依据，因此近年来望诊的客观化、现代化研究得到广大关注并取得了一定进展。但是，目前望诊客观化研究均是利用二维图像进行表面信息的采集，其中对舌象的表面形态以及舌表面的舌乳头大小分布等细节信息不能给出定量而形象的表示。因此在中医望诊的客观化研究领域，真实形象的表达其表面形态和细节信息有客观必要性。

发明内容

本方法针对望诊客观化研究中采用二维图像信息采集方式的局限性，提出基于光度立体技术进行望诊信息表面的三维重建。本发明将光度立体技术、图像摄取技术和数字图像处理技术等应用于中医望诊的客观化研究，实现中医望诊信息表面的三维重建。

光度立体技术是指通过光照变化对物体三维模型进行建模的方法，其基本原理是在至少三个不同光源方向的情形下，保持物体和拍摄仪器的位置不变，并且拍摄角度和焦距不变，在不同光源方向下拍摄多幅图像，利用图像的光强来计算物体表面的方向向量，并计算表面深度值进而得到物体表面的三维重建图。利用光度立体技术基于朗伯体反射模型计算物体表面的单位法向量，用以表达物体的细节几何，由单位法向量计算表面所有点的深度信息，并通过三维重建图形象的表达物体的表面形态和属性信息，从而达到物体三维表面的重建效果。

本发明具体技术方案如下：

基于光度立体技术的中医望诊信息表面三维重建方法，基于两组图像，这两组图像分别是以黑球为拍摄对象的高光黑球图像，以及以表面待重建对象为拍摄对象的表面待重建对象图像，其中每个光源下的高光黑球图像都与该光源下的表面待重建对象图像匹配，通过每一个高光黑球图像查找其高光点位置，获得光源方向，并作为配对的表面待重建对象图像的光源方向向量；然后利用光度立体技术计算得到待重建表面的单位法向量，再由代数法计算待重建表面所有点的深度信息，得到待重建对象的三维重建图。

上述两组图像的拍摄条件相同，拍摄条件为，在一间黑室内拍摄对象和拍摄仪器位置固定不变，并且拍摄角度和焦距不变，选取p个不同方向的光源依次照射在拍摄对象上，p大于等于3；所述的高光黑球图像即为上述拍摄条件下获取的p个不同光源方向下的高光黑球图像，其中每个图像均包括黑球部分和背景部分；所述的表面待重建对象图像即为上述拍摄条件下获取的p个不同光源方向下的表面待重建对象图像，其中每个图像均包括待重建部分和背景部分，获取两组图像后传至计算机。在计算机中完成图像的预处理；再利用光度立体技术实现待重建对象的表面三维重建，计算机处理具体过程包括以下步骤：

⑴对上述两组图像进行预处理，包括：

①图片校正：利用上述拍摄仪器在任意相同静态场景下重新采集一组不同曝光时间的图像，根据这组不同曝光时间的图像通过合成得到拍摄仪器的光电响应曲线，利用该光电响应曲线对在黑室内所拍摄的两组图像进行校正，得到校正后的两组图像；

②图片裁剪：对于校正后的两组图像，利用智能图像处理软件光影魔术手得到包含完整目标的矩形区域，分别把第一组高光黑球图像中包含完整黑球部分的同一区域裁剪出来，把第二组图像中包含完整待重建部分的同一区域裁剪出来，以减小图像中的背景区域，裁剪后每一组内各个图像大小相同、包含的区域相同；

③制作Mask图片：从裁剪后的两组图像中分别选出一张高光黑球图像和一张拍摄对象图像，利用智能图像处理软件Photoshop的磁性套索工具分别从高光黑球图像和拍摄对象图像中提取出黑球部分的轮廓和待重建部分的轮廓，然后将轮廓内的部分均用白色标记，其余部分均用黑色标记，从而得到两张新的图片，分别称为高光黑球的Mask图片和重建部分的Mask图片；

⑵计算光源的方向向量L：

①读取高光黑球的Mask图片并计算高光黑球的质心位置坐标和半径，然后依次读取裁剪后的第一组高光黑球图像并计算每个高光黑球图像的高光点位置坐标，所述的高光黑球的Mask图片与裁剪后的高光黑球图像大小相同，且Mask图片的质心位置坐标与裁剪后第一组高光黑球图像中各个图像的质心位置坐标相同；

②根据黑球高光点的镜面反射原理，计算归一化的光源方向向量L，计算公式如下：

L＝2(NR)N-R

其中，N＝(N₁,N₂,...,N_p)^T，N_i＝(x_i,y_i,r)，i为整数，且i∈[1,p]，分量x_i表示裁剪后的第i张高光黑球图像中高光点横坐标与高光黑球的质心横坐标之差；分量y_i表示裁剪后的第i张高光黑球图像中高光点纵坐标与高光黑球的质心纵坐标之差；分量r表示高光黑球的半径；R表示相机所在轴视角方向的单位向量；

⑶计算待重建对象的三维重建图中任意像素点q的单位法向量n，计算公式如下：

$n=\frac{1}{3}(n_R+n_G+n_B)$

其中，n_R,n_G,n_B分别表示待重建对象的三维重建图中q像素点的R、G、B三个通道对应的单位法向量，s分别代表R、G、B，g_s＝A_s ^-1b_s，A_s＝I_sL，L为归一化的光源方向向量，I_s＝(I_1qs，I_2qs，...，I_pqs)，I_iqs表示裁剪后第二组表面待重建图像中第i个图像的q像素点所对应的s值，其中i为整数且1≤i≤p；

⑷得到重建表面的法向量场，根据切平面理论：法向量与切平面的两个主轴垂直，用深度值表示表面切向量，对表面所有像素点的深度建立约束方程，计算待重建对象的三维重建图的深度信息矩阵Z，所述的三维重建图与裁剪后的第二组图像大小相等，均为U×V，深度信息Z的计算公式为：

Z＝C^-1d，

其中，C表示一个(2(U×V-1),U×V)的二维稀疏矩阵，具体为

其中，c(i,j)表示矩阵C中第i行第j列的元素，i、j为整数，满足1≤i≤2(U×V-1),1≤j≤U×V；

d表示(2(U×V-1),1)的矩阵，表示如下：

其中，d(k,1)表示矩阵d中第k行第1列的元素；

n_X表示(U×V-1,1)的矩阵，具体计算方法如下：首先将步骤（3）得到的待重建对象的三维重建图中最后一行像素点法向量在X轴上的分量n_x均设为0，然后以列递增的顺序依次提取待重建对象三维重建图中所有像素点单位法向量的n_x分量，并组成列向量，提取前U×V-1个元素组成一个(U×V-1,1)的矩阵，即为n_X；

n_Y表示(U×V-1,1)的矩阵，具体计算方法如下：首先将步骤（3）得到的待重建对象的三维重建图中最后一列像素点法向量在Y轴上的分量n_y均设为0，然后以列递增的顺序依次提取待重建对象三维重建图中所有像素点单位法向量的n_y分量，并组成列向量，提取前U×V-1个元素组成一个(U×V-1,1)的矩阵，即为n_Y；

⑸建立深度信息Z与待重建对象的三维重建图中各像素点的对应关系，绘制待重建对象的三维重建图；所述的深度信息Z与待重建对象的三维重建图中的像素点的对应关系为：从包含U×V个元素的深度信息Z中依次提取U个元素作为一列，顺序组成U行V列的矩阵，该矩阵中的深度信息分别与大小为U×V的待重建对象的三维重建图中的像素点一一对应。

最后利用待重建表面的深度信息，绘制待重建表面的三维重建图，实现三维表面的重建。

有益效果

本方法可实现望诊信息表面的三维重建，即从重建结果中形象的观察到舌象或面部的表面形态和细节属性特征，并利用重建的三维视图多角度观察其表面形态。本方法首先计算待重建表面的单位法向量；然后再计算待重建表面每一点的深度信息，最后输出待重建表面的三维重建图。从对整个舌象的重建效果可以看出舌体的整体形态，可以定量给出舌体上的裂纹等信息；并且在三维重建图中可以从各个角度进行观察，更形象的表达真实舌体的完整信息。利用舌象表面的三维重建图进行中医舌象的信息采集和表达，使其表达信息更完整形象，有利于推动中医舌诊的客观化研究进程。

附图说明

图1：整体方法流程图；

图2：由光度立体技术求解法向量流程图；

图3：代数法求解深度信息流程图；

具体实施方式

本实施例的拍摄条件是在一间黑室内，背景用黑绒布遮掩，数码相机放在舌体模型正前方，数码相机和舌体模型的位置固定不变并且拍摄角度和焦距不变，取8个LED光源沿直径方向均匀分布于半径为28厘米的竖直圆周上，数码相机位于其圆周中心并与舌体模型等高，舌体模型位于数码相机正前方20厘米处。依次点亮8个LED光源，通过软件操作数码相机分别在不同光源方向下拍摄一张表面待重建舌象图像，共计8张。然后在同样的拍摄条件下，将舌体模型换成黑球，拍摄8张高光黑球图像，然后将这两组图像输入至计算机进行预处理，图像格式为JPGE图像，并利用光度立体技术进行重建，整个方法流程图见图1，具体实施过程如下：

1.对拍摄的表面待重建舌象图像和高光黑球图像进行预处理，具体过程包括如下步骤：

⑴图片校正：利用上述拍摄仪器拍摄一组不同曝光时间的图像，曝光时间依次为H=[2,1,0.5,0.25,0.125,0.0667,0.0333,0.0167]，根据公式g(I_mt)＝lnE_m+lnH_t拟合得到该数码相机的光电响应曲线g，I_mt表示第t个曝光时间图像的m位置处像素点的像素值，1≤t≤8，t为整数，其中g表示相机的光电响应曲线，即像素值与光亮值之间的对应关系；E_m表示辐照度值，H_t表示第t个曝光时间，然后根据该数码相机的光电响应曲线g，将上述拍摄到的表面待重建舌象图像和高光黑球图像中每个图像的像素值映射成光亮值，理想的光电响应曲线即像素值与光亮值是线性关系，根据理想的线性关系再将光亮值按比例映射回0-255的像素值，这样即可得到校正后的像素值；对于R、G、B三个颜色通道各自的响应曲线，使用相同的方法，对所有图像进行校正，即可得到校正后的两组图像；

⑵图片裁剪：对于校正后的两组图像，利用智能图像处理软件光影魔术手分别把第一组高光黑球图像中包含完整黑球部分的同一区域裁剪出来；再利用光影魔术手把第二组图像中包含完整重建舌象的同一区域裁剪出来，以减小图像中的背景区域，裁剪后每一组内各个图像大小相同、包含的区域相同，且保证裁剪后每个光源下两组图像的一一对应关系；

⑶制作Mask图片：从裁剪后的两组图像中选出一张高光黑球图像和一张拍摄舌象图像，用Photoshop的磁性套索工具分别从这两张高光黑球图像和拍摄舌象图像中提取出高光黑球部分的轮廓和重建舌象的轮廓，然后将轮廓内的部分均用白色标记，其余部分均用黑色标记，从而得到两张新的图片，分别称为高光黑球的Mask图片和重建舌象的Mask图片；

2.求解各个光源的方向向量，具体步骤如下：

⑴读取高光黑球的Mask图片，先对图像进行灰度化，然后计算图片中高光黑球的质心位置坐标和半径，其中图像上点的坐标就是该点所在的行和列的索引值，该高光黑球的Mask图片与裁剪后第一组高光黑球图像大小相同，因此该高光黑球Mask图片中高光黑球的质心坐标与裁剪后第一组高光黑球图像中各个图像的质心位置坐标相同；

⑵依次读取8个不同光源方向的高光黑球图像，设置R、G、B的阈值为220，对于每个高光黑球图像，得到R、G、B大于阈值的高光区域，再利用求质心的方法计算高光区域的质心，即作为各个高光点的位置坐标；

⑶计算N，N＝(N₁,N₂,...,N_p)^T，N_i＝(x_i,y_i,r)，i、p为正整数，p为权利要求1中所述的不同方向光源的个数，且i∈[1,p]，分量x_i表示裁剪后的第i张高光黑球图像中高光点横坐标与高光黑球的质心横坐标之差；分量y_i表示裁剪后的第i张高光黑球图像中高光点纵坐标与高光黑球的质心纵坐标之差；分量r表示高光黑球的半径；

⑷根据公式L＝2(NR)N-R，计算得到归一化的光源方向向量L，其中，相机所在轴视角方向的单位向量R＝(0,0,-1)；

3.利用光度立体技术计算待重建对象的三维重建图中任意像素点q的单位法向量，流程图见图2，具体过程如下：

(1)依次读取8个预处理后的待重建舌象图像，得到所有图像个点的像素值，取8个待重建舌象图像中像素点q的R、G、B值，分别组成I_s，其中s分别代表R、G、B，I_s＝(I_1qs，I_2qs，...，I_pqs)，I_iqs表示裁剪后第二组图像中第i个图像的任意像素点q所对应的s值，其中i为整数且1≤i≤p，p的定义如前文所述；

(2)由像素值I_s计算A_s、b_s：A_s＝I_sL，其中L为归一化的光源方向向量；

(3)由A_s,b_s根据公式g_s＝A_s ^-1b_s计算得到g_s；

(4)计算待重建舌象的三维重建图中像素点q的R、G、B三个通道分别对应的单位法向量，计算公式为：其中s分别代表R、G、B；

(5)得到待重建舌象三维重建图中任意像素点q的R、G、B三个通道对应的单位法向量n_R,n_G,n_B，取其平均值作为该像素点q的单位法向量n，即待重建舌象的三维重建图中任意像素点q的单位法向量n的计算公式如下：

$n=\frac{1}{3}(n_R+n_G+n_B)$

逐像素点计算即可得到待重建舌象三维重建图中所有像素点的单位法向量；

4.计算待重建舌象的三维重建图中各像素点的深度信息Z，相应的流程图见图3，所述的三维重建图与裁剪后的第二组图像大小相等，均为U×V，步骤如下：

⑴建立一个(2(U×V-1),U×V)的二维稀疏矩阵，即矩阵C，具体为：

其中，c(i,j)表示矩阵C中第i行第j列的元素，i、j为整数，满足

1≤i≤2(U×V-1),1≤j≤U×V；

举例说明，当U＝2,V＝3时，矩阵C具体形式如下：

$C_{\mathrm{10,6}}=\left[\begin{array}{cccccc}1& -1& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 1& -1& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 1& -1\\ 1& 0& -1& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 1& 0& -1& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right]$

⑵计算n_X以及n_Y：

n_X表示(U×V-1,1)的矩阵，具体计算方法如下：首先将步骤3得到的待重建舌象的三维重建图中最后一行像素点法向量在X轴上的分量n_x均设为0，然后以列递增的顺序依次提取待重建舌象三维重建图中所有像素点单位法向量的n_x分量，并组成列向量，提取前U×V-1个元素组成一个(U×V-1,1)的矩阵，即为n_X，

假设三维重建图大小为U×V，U＝2,V＝3，待重建舌象三维重建图中所有像素点单位法向量在X轴上的分量 $n_x=\left[\begin{array}{ccc}1& 0& 0\\ 1& 1& 1\end{array}\right],$ 令三维重建图中最后一行像素点法向量在X轴上的分量n_x均设为0，即 $n_x=\left[\begin{array}{ccc}1& 0& 0\\ 0& 0& 0\end{array}\right],$ 以列递增的顺序提取待重建舌象三维重建图中所有像素点单位法向量的n_x分量，并组成列向量即 $\left[\begin{array}{c}1\\ 0\\ 0\\ 0\\ 0\\ 0\end{array}\right],$ 提取前U×V-1个元素组成一个(U×V-1,1)的矩阵，即为 $n_X=\left[\begin{array}{c}1\\ 0\\ 0\\ 0\\ 0\end{array}\right]$

n_Y表示(U×V-1,1)的矩阵，具体计算方法如下：首先将步骤（3）得到的待重建舌象的三维重建图中最后一列像素点法向量在Y轴上的分量n_y均设为0，然后以列递增的顺序依次提取待重建舌象三维重建图中所有像素点单位法向量的n_y分量，并组成列向量，提取前U×V-1个元素组成一个(U×V-1,1)的矩阵，即为n_Y；

⑶根据n_X以及n_Y，建立一个(2(U×V-1),1)的d矩阵，表示如下：

其中，d(k,1)表示矩阵d中第k行第1列的元素

⑷由矩阵C和矩阵d计算待重建舌象三维重建图的深度信息矩阵Z，其计算公式如下：

Z＝C^-1d

然后，建立深度信息矩阵Z与待重建舌象的三维重建图中各像素点的对应关系，即从包含U×V个元素的深度信息矩阵Z中按行递增的顺序依次提取U个元素作为一列，顺序组成U行V列的矩阵，该矩阵中的深度信息分别与大小为U×V的待重建舌象的三维重建图中的像素点一一对应，然后绘制待重建舌象的三维重建图。

例如，三维重建图大小为U×V，U＝2,V＝3，计算得到深度信息矩阵 $Z=\left[\begin{array}{c}1\\ 1\\ 0\\ 1\\ 1\\ 0\end{array}\right],$ 按行递增的顺序依次提取2个元素作为一列，顺序组成2行3列的矩阵 $\left[\begin{array}{ccc}1& 0& 1\\ 1& 1& 0\end{array}\right],$ 该矩阵中的深度信息分别与大小为2×3的待重建舌象的三维重建图中的像素点一一对应，然后根据所有像素点所在行和列以及深度信息，利用MATLAB的三维绘图即可绘制待重建舌象的三维重建图。

Claims (2)

1.基于光度立体技术的中医望诊信息表面三维重建方法，基于两组图像，这两组图像分别是以黑球为拍摄对象的高光黑球图像，以及以表面待重建对象为拍摄对象的表面待重建对象图像，两组图像的拍摄条件相同，选取p个不同方向的光源依次照射在拍摄对象上，p大于等于3；所述的高光黑球图像即为上述拍摄条件下获取的p个不同光源方向下的高光黑球图像，其中每个图像均包括黑球部分和背景部分；所述的表面待重建对象图像即为上述拍摄条件下获取的p个不同光源方向下的表面待重建对象图像，其中每个图像均包括待重建部分和背景部分，两组图像获取后传至计算机，其特征在于还包括以下步骤：

⑴对上述两组图像进行预处理，包括：

①图片校正：根据拍摄仪器的光电响应曲线对上述高光黑球图像和表面待重建对象图像进行校正，得到校正后的两组图像；

②图片裁剪：对于校正后的两组图像，分别把第一组高光黑球图像中包含完整黑球部分的同一区域裁剪出来，把第二组图像中包含完整待重建部分的同一区域裁剪出来，以减小图像中的背景区域，裁剪后每一组内各个图像大小相同、包含的区域相同；

③制作Mask图片：从裁剪后的两组图像中选出一张高光黑球图像和一张拍摄对象图像，从高光黑球图像中提取黑球部分的轮廓，从拍摄对象图像中提取待重建部分的轮廓，然后将轮廓内的部分均用白色标记，其余部分均用黑色标记，从而得到两张新的图片，分别称为高光黑球的Mask图片和待重建部分的Mask图片；

⑵计算光源的方向向量L：

①读取高光黑球的Mask图片并得到高光黑球的质心位置坐标和半径，然后依次读取裁剪后的高光黑球图像并得到每个高光黑球图像的高光点位置坐标；

②根据黑球高光点的镜面反射原理，计算归一化的光源方向向量L，计算公式如下：

L＝2(NR)N-R

其中，N＝(N₁,N₂,...,N_p)^T，N_i＝(x_i,y_i,r)，i、p为正整数，p为上述的不同方向光源的个数，且i∈[1,p]，分量x_i表示裁剪后的第i张高光黑球图像中高光点横坐标与高光黑球的质心横坐标之差；分量y_i表示裁剪后的第i张高光黑球图像中高光点纵坐标与高光黑球的质心纵坐标之差；分量r表示高光黑球的半径；R表示相机所在轴视角方向的单位向量；

⑶分别计算待重建对象的三维重建图中任意像素点q的单位法向量n，计算公式如下：

$n=\frac{1}{3}(n_R+n_G+n_B)$

其中：n_R,n_G,n_B分别表示待重建对象的三维重建图中像素点q的R、G、B三个通道对应的单位法向量，s分别代表R、G、B，g_s＝A_s ^-1b_s，A_s＝I_sL，L为归一化的光源方向向量，I_s＝(I_1qs,I_2qs，...，I_pqs)，I_iqs表示裁剪后第二组图像中第i个图像的q像素点所对应的s值，其中i为整数且1≤i≤p；

⑷计算待重建对象的三维重建图的深度信息矩阵Z，所述的三维重建图与裁剪后的第二组图像大小相等，均为U×V，深度信息矩阵Z的计算公式为：

Z＝C^-1d，

其中，C表示一个(2(U×V-1),U×V)的二维稀疏矩阵，具体为

其中，c(i,j)表示矩阵C中第i行第j列的元素，i、j为整数，满足1≤i≤2(U×V-1),1≤j≤U×V；

d表示(2(U×V-1),1)的矩阵，表示如下：

其中，d(k,1)表示矩阵d中第k行第1列的元素，

n_X表示(U×V-1,1)的矩阵，具体计算方法如下：首先将步骤（3）得到的待重建对象的三维重建图中最后一行像素点法向量在X轴上的分量n_x均设为0，然后以列递增的顺序提取待重建对象三维重建图中所有像素点单位法向量的n_x分量，并组成列向量，提取前U×V-1个元素组成一个(U×V-1,1)的矩阵，即为n_X，

n_Y表示(U×V-1,1)的矩阵，具体计算方法如下：首先将步骤（3）得到的待重建对象的三维重建图中最后一列像素点法向量在Y轴上的分量n_y均设为0，然后以列递增的顺序提取待重建对象三维重建图中所有像素点单位法向量的n_y分量，并组成列向量，提取前U×V-1个元素组成一个(U×V-1,1)的矩阵，即为n_Y；

⑸建立深度信息矩阵Z与待重建对象的三维重建图中各像素点的对应关系，绘制待重建对象的三维重建图；所述的深度信息矩阵Z与待重建对象的三维重建图中像素点的对应关系为：从包含U×V个元素的深度信息矩阵Z中依次提取U个元素作为一列，顺序组成U行V列的矩阵，该矩阵中的深度信息分别与大小为U×V的待重建对象的三维重建图中的像素点一一对应。

2.根据权利要求1所述的基于光度立体技术的中医望诊信息表面三维重建方法，其特征在于：所述的拍摄条件为，在一间黑室内拍摄对象和拍摄仪器位置固定不变，并且拍摄角度和焦距不变。