专利名称：基于物理的植物叶子模拟方法
技术领域：
本发明属于计算机真实感图形学、虚拟现实、几何造型、计算机动画、计算机三维建模和植物仿真学相结合的领域，具体涉及一种通过物理模拟对植物叶子造型的方法。
背景技术：
在计算机图形学中，对于自然场景中植物的建模和绘制一直是绘制技术发展过程中的重点。对植物进行造型、材质表达以及真实感绘制，一直是极富挑战性的研究工作。
众所周知，植物往往由很多具有不同外形特征的器官组成，如茎、叶子、花、果实等。而自然界中的植物，与大部分人造形体的最大区别在于人造形体往往按照人们的审美和应用需求，在制造过程中，其几何形状是可以控制，并易于被描述和掌握；而对于自然景物的植物而言，其几何形状的形成，则是自然界中各种复杂内外因素，如外部病虫害入侵、气候因素、人为影响等，经过长时间互相影响、共同作用的结果。因此，作为自然界中的植物，通常拥有非常特殊的形态。尤其植物的叶子，随着季节与气候的变化，其形态变化十分显著，这也在相当程度上决定了植物外观特征。而在以往的计算机辅助建模中，用简单、固定的几何图形很难真实地表现叶子在不同环境和季节下的外形特征。同时，由于现有技术的限制，在计算机图形学真实感绘制领域中，大部分绘制出来的场景均为静态场景，即场景中的物体几何形状并不能随时间的变化而变化，这不仅与自然界的真实情况不符，绘制的真实感不强，也极大地限制了植物仿真技术的进一步发展。因此，需要针对作为自然场景中最主要组成部分之一的植物树叶的几何形态的特征，以及它们随时间变化而变化的情况进行了分析，使叶子的绘制结果更加贴近自然真实。通过这些的研究，就能在虚拟现实的相关应用中构建更加逼真的场景，也可以让植物的真实感绘制在游戏、影视、生物、环保等领域得到进一步应用。
对于自然物体外观仿真在近四十年来自然场景中植物的模拟和绘制得到了广泛的应用。此前关于植物建模方面的研究工作主要包括Honda[1]、Webber[2]等提出的参数化建模方法；Reeves[3]提出的基于粒子系统的仿真；Bloomenthal[4]、Holton[5]基于spline骨架的造型；Fournier[6]、Viennot等[7]提出的基于随机系统的建模方法；Lintermann和Deussen[8]提出的基于图结构的建模；De Reffye等[9]、Prusinkiewicz和Lindenmayer[10]提出的重写(rewriting)规则的系统；Lefebvre和Nyeret[11]、Federl和Prusinkiewicz[12]提出的基于物理的植物纹理模拟。其中Lintermann和Deussen在他们的Xfrog交互式植物模型系统中结合了植物的位置信息，沿轴向方向采样点的位置信息可以影响节点和枝条的长度、分支的角度和数量，以及它们的特征。
Sung[16]提出过一种通过三角形网格构建植物叶子表面，并通过用户交互界面来编辑植物叶子变形的方法，其生成的模型结果具有很强的真实感，但是其方法缺乏很好的植物学和物理学依据，变形结果由人为设定。
直接对于叶子形态进行构型的工作目前还不多见，Judd[14]等人给出了一个描述叶子形状的表达方法，此后被广泛应用。而对于植物叶脉的造型过程，Dengler和Kang[13]的研究表明，叶脉的模型与叶子的形状具有强烈的相关性，因此必须考虑叶脉与树叶形状之间的关系。2005年的SIGGRAPH会议上，Runions和Prusinkiewicz等人[15]提出了一种更加有效的叶脉构造方法，用以生成与叶子外形相关的叶脉模型。基于对生长中叶子的植物学特性观察，这种方法在叶肉模型表面上均匀地布置了一些营养点。然后，假设叶脉的生长的目的是到达每一个营养点，并由此可为叶肉上的每个细胞输送水分和营养，进而决定整个叶脉的生长的方式。用这种叶脉模型生成方法可以得到真实感很强的叶脉构型，不过Runions方法的不足之处在于，其生成的叶脉结构仅仅在空间分布和生长形态上与真实叶脉接近，但其生成的结果仅是叶脉的平面效果，即是一个二维纹理。
因此，在计算机图形学的很多应用场合，尤其是真实感绘制领域，急需一种高效、逼真的方法来生成具有高度真实感的植物叶子几何模型。
发明内容
本发明的目的在于运用计算机真实感图形学绘制和虚拟现实技术，提供一种对植物叶子进行模拟的方法，即植物叶子不同季节和时期动态建模和绘制的方法。
如前所述，目前的一些植物建模大部分是基于语义和分形的构型和造型技术，它们的共同特点就是实现简单，生成速度快。另外一些基于图像的植物建模技术可以生成相当拟真的植物几何模型。但是，这些现有方法的共同问题就是不能很好模拟植物的动态效果，难以表现出植物在不同环境下的形态变化，如风吹动时候的摇曳效果，秋季到来植物叶子的老化卷曲等变化效果。本发明通过对植物细胞力学结构的研究，确立了构成植物的最基本单位细胞是由较硬的细胞壁和较软的内部结构组成的分析结果，并通过这些对于植物组织组成的分析，构建了植物的力学模型对植物的形态进行模拟，从而对植物在不同时态和季节下的形态和变化，进行了很好的模拟。
本发明对植物叶子进行模拟的方法主要分为5个步骤植物叶子轮廓提取、植物叶脉构造、叶子表面模型三角化、构造叶子质点一弹簧双层模型结构、进行物理模拟运算。下面对本发明的技术方案进行详细描述 1.植物叶子轮廓的提取 构造植物叶子的三维几何模型，首先要从真实世界中的植物模型进行采样，获得一些实际的数据，以便对植物叶子的模拟更加接近真实世界中的情况。对植物叶子轮廓的提取过程具体包括 1.1叶子图像二值化 首先对植物叶子进行图像扫描，得到叶子的平面图像。然后对图像进行二值化，也就是将扫描后的叶子平面图像上所有的像素点的取值定为黑色或者白色，这样有利于区分哪些像素属于叶子内部，哪些像素不属于。图像二值化之后，叶子平面图像变为一副黑白二值图像。二值图像可以显著减少在获取叶子轮廓时候的运算量。通过设置不同的分界线取值，可以得到不同的二值图像，如图1所示。
1.2叶子轮廓的步进方格提取法(Marching square) 本方法所使用的步进方格法(Marching Square)是著名的步进立方体算法(MarchingCube)[19][20]在二维情况下的一种修改过后的方法。步进方格法以叶子轮廓可以由多边形表示为基本假定，同时也假定这个多边形也是一个闭合曲线。用方格对叶子二值图像从左上到右下的顺序进行扫描线方式的水平扫描，直到方格的某个角进入叶子内部的区域，即方格经过叶子边界，这时，方格停止扫描线方式的水平扫描，转而开始沿着叶子的边界以逆时针方向运动。由于叶子的轮廓的封闭和连续假设，当方格的大小合适时，能够得到非常逼近现实中叶子轮廓的多边形轮廓。
在通过扫描叶子获得了叶子的二值图像之后，就可以开始使用步进方格法(MarchingSquare)获取叶子的多边形轮廓了。在步进方格法中，一旦方格的四个角的像素中有黑色的像素，则说明方格开始横跨叶子内部和外部之间的区域。通过方格四角像素的颜色，可以判别方格的哪条边被叶子轮廓所经过，如图2所示，就可以得到方格中叶子边的走向和位置，从而得到多边形叶子轮廓的一小段，步进方格法就是将获得的轮廓线段的每一小段找到并连接起来的过程。通过步进方格对整个叶子图像的扫描，就可以得到叶子的多边形轮廓。提取方式示意图如图2所示，详细过程可以参考文献[16]。
在二值图的步进方格法(Marching Square)中，由于每个像素有黑白两种取值，方格四角颜色可能出现的情况有24共16种，如图3所示，通过这16种方格的情况，就可以判定方格位于叶子轮廓的位置，从而对轮廓进行提取。不同大小的方格能够提取不同细节的多边形叶子轮廓，如图4所示，图的上部是一片枫叶的背面图像，图的下部是枫叶方框中的部分通过步进方格法提取轮廓，方格边长选择不同大小时侯得到的结果，从左到右方格所选边长依次是8像素、5像素、3像素、2像素，可以看出，方格越小所提取的轮廓越精细。
2.植物叶脉构造 在真实的自然环境中，植物叶子的变形会受到多种多样的内外因素的影响，为了方便研究且不失一般性，本发明考虑其中最主要的影响因素——叶子的水分供给。因为在叶子的老化、枯萎的过程中，叶脉的模型对叶子的几何构型有着非常重要的影响。而建立逼真的叶脉分布是真实模拟植物叶子变形的基础。
对此，Runions[15]曾提出了一种有效的生成叶脉的方法，一般被称为“基于营养点的叶脉生成方法”。该方法根据于叶子生长中所表现出来的一些植物学现象，在叶子的表面均匀地分布一些代表营养物质的“营养点(或称为生长素点)”，并由此假设叶脉的生长的目标就是到达并连接叶子上的各个营养点，从而为生长中的叶肉上的每个细胞运送水分和营养物质。这样，叶片的外形和营养点的分布就决定了叶脉生长的方式，这样就很好地把叶脉的构型和叶子形状有机联系起来。
具体步骤如图5所示(图5中实心圆点表示营养点，环状点表示叶脉上的节点) 1.在已有的叶子轮廓的内部随机散布营养点； 2.为每个营养点寻找最近的叶脉节点； 3.寻找成为营养点最近叶脉节点的这些节点，并记录之； 4.为每个这样的“最近节点”计算生长方向，方向由对其产生作用的营养点对其影响的加权平均和计算出； 5.按照计算出的方向以限定间距生长出新的叶脉节点； 6.根据给定的某个半径，计算营养点是否被叶脉达到； 7.删除已经达到的营养点； 8.轮廓调整； 9.随机散布新的营养点； 10.删除距离已有营养点过近、以及距离叶脉节点过近的新营养点； 11.回到步骤1，重新计算下一个过程中叶脉应该如何生长。
经过多轮循环处理，所生成的叶脉实际上是由连接各叶脉节点的线段组成。这种现有的技术可以用来真实地生成很多种类叶子所具有的叶脉结构，包括玉兰叶、樱桃叶、枫叶和梧桐叶等形状的叶子。虽然L-system[10]也可以生成类似的叶脉结构，但往往缺乏植物学上的理论依据。
在植物学中，绝大多数双子叶植物的叶子，在生长的初期，叶子中连接叶根部和叶顶端的中脉就已经形成，而肋状的分支叶脉都是从中脉分化生长出来，因此中脉的物理强度和几何直径都远大于其他分支叶脉。同样，在趾状叶子生长过程中，每个趾结构会形成一个自身的中脉，并逐渐生长出基于该中脉的分支叶脉结构。利用这种影响，本发明改进了现有的叶脉生成方法，具体内容即在叶脉生成的过程中，提高模拟过程中主脉和中脉的强度和生长速度，因为真实叶子的主叶脉能够达到和叶子轮廓基本相同的生长速度，而比次级叶脉的生长速度快，一般设定主脉和中脉的生长速度是次级叶脉的1.5-4倍，优选为2倍。这样，就可以确保所生成的叶脉结构能够有更现实的叶脉结构物理强度。此外，生长速度很快的中脉和主脉可以通过它们的空间分布影响次级叶脉的生长方式。通过调整和测试，本发明可以获得叶脉结构和形状分布更加逼真的效果。如图6所示，最左侧为使用原方法生成的叶脉，中间为本发明的方法生成的叶脉，最右侧为真实的叶子图像，可以看出本方法生成的叶脉更加接近真实情况。
在后面的物理模拟过程中，被抽象成为线段的每段叶脉的硬度直接影响着物理模拟变形的结果。而叶脉的强度主要有叶脉的粗细(即直径)决定，本发明根据Murray[17]提出的Murray定律来对结构中叶脉的直径进行设置 这里d表示叶脉段的直径，也就是叶脉的宽度，dparent、dchild1和dchild2分别表示叶脉分枝处分枝前的叶脉段的直径和分枝后形成的两个叶脉段的直径，n为常数。Runions[15]在他的工作中用n＝3来求得叶脉宽度进行叶子的绘制。而在本发明中，叶脉段直径的设置并不是为了在本步骤中绘制出一段宽度的叶脉段，而是为后续的物理模拟过程中叶脉强度的设定做准备，由于在物理模拟的过程中，假定每段叶脉的直径d是一个不变的常量，最初的叶脉硬度成比例于叶脉的横截面面积π(d/2)2，因此，经过实验发现，在本发明中n＝1或者n＝2更加准确一些，更加接近真实。这样，每段叶脉最终的叶脉直径d将会直接体现为叶脉的硬度，而参与到后面的物理模拟运算中去。
3.叶子表面模型三角化 仅仅拥有叶子的叶脉和轮廓，还不足以构建叶子的物理结构，对叶子进行物理模拟。因此需要构造整个叶子的三维物理模型。
作为植物的一种器官，叶子的表面是由很多细胞组成的。并且分为上表皮、栅栏组织、海绵组织、下表皮等竖直方向上的结构。但是在对叶子进行宏观的物理模拟的时候，由于叶子的厚度相对显得非常的小，要远远小于叶子在水平方向上的尺度，因此可以将叶子表面的的结构在宏观上抽象简化成为一个二维平面。而为了能够更好的表达这个二维平面上叶子的各个部分之间的连接关系和相互作用，本发明使用三角形网格来构建叶子叶肉表面的结构。
本发明使用的三角化方法是Delaunay三角化方法[21][22][23]。由于已经拥有叶脉和叶子轮廓作为三角形网格的基础和限制，这里选择受限制下的Delaunay三角化方法(Constrained Delauany Triangulation)来对本发明的模型进行三角网格化，即以已有的叶脉线段和轮廓线段为约束，进行三角化的过程。在对叶子表面进行三角化的过程中，本发明使用了CGAL开源程序进行辅助，并通过限制每个三角形的最小角的角度和最小面积来使得三角化后的叶子表面三角网格更加均匀，从而利于下面步骤中细胞力学模型的构建，能够很好的支持接下来的物理模拟过程。三角化过程中的具体步骤是使用CGAL库提供的Exact_predicates_inexact_constructions_ kernel类、Delaunay_mesher_2类、Delaunay_mesh_face_base_2类、Constrained_Delaunay_triangulation_2类这4种数据结构来组织以叶子轮廓和叶脉作为限制的线段和点数据，并使用refine_Delaunay_mesh_2函数进行受限制的Delaunay三角化，数据结构Criteria进行三角形最小角度和最小面积控制，一般选取参数Criteria(0.005，0.5)，这时候表示三角形最小边的边长不能小于0.005，最小角的角度不能小于0.5弧度，从而使得生成的三角形网格更加适于物理模拟的过程。部分三角化结果如图7所示。
4.构造叶子质点-弹簧双层模型结构 在拥有植物叶子的轮廓和叶脉之后，就可以考虑组织构建叶子的物理结构了。
植物的树叶在不同的季节展现出来的是不同的形态，当夏季阳光和水分非常充足的时候，树叶会伸展开，呈现上表面凸出的形态。而秋季来到时，很多阔叶植物的叶子会由于气温的下降和水分的流失而向内卷曲，空气中湿度的变化、植物体内水分的变化以及其他一些条件的变化都会使叶子在几何形状上发生有规律的变化。
植物的细胞会由于外界环境的变化而变形，而植物叶子的变形正是叶子上众多细胞共同变形的结果。通过对植物叶子的观察，本发明发现，在微观尺度上，绝大部分植物的叶子主要由上表面的叶肉和大部分位于叶子下表面的树状结构的叶脉组成。叶肉主要由上表皮、栅栏组织、海绵组织和下表皮组成，叶脉主要由众多维管束组成。
维管束和叶肉细胞有着十分不同的硬度和物理结构。叶肉细胞在水分流失的时候会有很明显的体积收缩，而在相同情况下，维管束的体积变化却相对小得多。这主要是由于维管束的细胞之间含有很多的木质和纤维。
当季节发生变化的时候，如秋季，随着植物叶水分的散失，物理结构不同的叶脉和叶肉就会产生不同程度的变形由于叶脉具有很多沿叶脉径向的维管束，这些径向的维管束中含有大量径向的纤维结构，就会导致在水分流失过程中，叶脉的收缩比率并不大；而叶肉由大量的饱含水分的叶肉细胞组成，细胞与细胞之间并没有强有力的纤维结构来保持形态不变，因此随着水分在秋季的流失，叶肉部分的细胞在沿叶肉平面的切向方向面积大幅收缩。
这种不相等的收缩比例在叶子的每个微小结构中都发生着，包括叶脉到达的每个位置。最终它会导致叶子这个局部向叶子上表面方向弯曲自身，即三维空间中几何形态的变化。最终变形的结果，也是叶子所有部位变形受力后共同作用的结果。
4.1细胞力学模型 本发明提出采用弹性力学模型来简化模拟叶肉细胞的物理受力情况。假定叶子上的一个细胞主要受到两种力，一种是由细胞质所产生的细胞内部的膨胀力fin，另外一种是由细胞壁产生的细胞外部的张力fout，如图8所示，其中最左侧的A图是对单个细胞受力的分析，中间的B图是多个细胞受力分析和弹簧模拟情况，最右侧的C图是抽象出来的质点弹簧网络模拟细胞物理结构的示意图。
为了更精确地描述叶子细胞所受的这两种物理力，构建收缩弹簧来表示细胞的膨胀力，同理，采用拉伸弹簧来表示细胞之间的张力。由此可知，一个细胞处于稳定几何形态的条件是这两种力相互平衡，可用下面方程表达 fin＝-fout＝-k(s-s0)(1) 其中，s和s0分别表示细胞饱满状态下的弹簧静止长度和细胞在失水收缩之后的弹簧长度，而收缩之后的长度s0也恰恰决定细胞的大小；k表示的是弹簧的弹性系数。通过实际观察发现，细胞内的液泡会随着秋季得到来，由于水分流失而显著缩小，由此就导致了细胞内部膨胀力fin的减少，细胞弹簧的静止长度s也随之减小。也就是说，随着细胞水分流失这个过程，细胞的内部组成的变化导致了其物理学和力学特性的变化，因此可以将其看作模拟细胞的弹簧的内在力学性质和参数发生变化的过程，而这个过程也可以看作弹簧变化的过程。叶肉细胞失去水分之后，其收缩的过程，在本发明的模拟方法中就表现为弹簧静止长度s减小的过程。
按照上面的简化模型，叶肉细胞组织可以由图8中C图所示的质点-弹簧网状结构来模拟。基于如上的植物学观察，本发明提出一种基于物理模拟的双层力学模型DLM(doublelayered model)来模拟叶子的变形。构建双层模型的基本思想源于一个自然界里的基本事实绝大部分叶脉是位于叶子背面，尤其对于一些阔叶植物来说，并且叶脉也凸出于叶子表面。而在双层模型中，本发明采用上层网格表示收缩率很高的叶肉上的细胞，而下层结构则用来表示非常硬且不易收缩的叶脉上的细胞。同时这样的两层结构模型还通过互相连接来保证力学上的稳定性，连接的弹簧结构代表了叶肉和叶脉细胞相连的结构。
4.2上层模型 为了构造模型中上层的质点-弹簧网络，由上述步骤中获得叶子轮廓，并对叶子表面生成叶脉，然后以叶脉和轮廓作为约束点进行Delaunay三角化，从而得到一张三角网格结构，网格的顶点表示细胞，作为质点，而网格中的边表示连接相邻质点的弹簧，用以表达细胞的力学结构。由此可以构建这样一个质点-弹簧网络，其中的弹簧是保证叶片形状的基础结构。
对于模型中，任何两个相邻的，并由弹簧相连接的质点i和质点j，用sij来表示连接质点i和质点j的一段弹簧，其中pi和pj相应表示质点i和质点j这两个质点的位置。通过对于质点i受力分析，弹簧sij在其上的作用fij可以由胡克定律得到 其中lij表示从pi到pj的向量，hij表示弹簧sij的弹性系数，dij表示弹簧sij的阻尼系数，rij表示弹簧sij的静止长度。这里的静止长度rij在模型中是一组随时间变化的参数。
表示质点i和j的相对运动速度，即lij对时间的导数。这里还可以知道vi和vj分别表示质点i和j的速度。
在方程(2)中，弹簧力的大小成比例于两质点实际长度与静止长度的差，阻尼力成比例于两个质点i和j的相对速度。这样，一个弹簧所连接的两个质点，受到大小相同，方向相反的力。弹簧的阻尼作用对称地作用于弹簧所连接的两个质点上，并不影响两个质点共同质心的运动。
质点i与相邻的若干质点j相连，作用在质点i上的力是由这些分别与i通过弹簧相连接的质点们的合力决定，这里Si表示的是所有与质点i相连的弹簧的集合，pj表示质点j的位置。这样，本发明建立了一个静止长度rij由方程(1)中变量s决定，并且弹性系数hij随时间t变化的质点-弹簧网络系统。最终上层模型的模拟方程可以描述为 方程(3)中，j表示与i相邻的并通过弹簧与i相连接的质点，sij表示连接质点i和质点j的一段弹簧，pi和pj相应表示质点i和质点j的位置，Si表示所有与质点i相连的弹簧的集合，lij表示从pi到位置pj的向量，hij(t)表示弹簧sij的弹性系数，dij表示弹簧sij的阻尼系数，rij(t)表示弹簧sij的静止长度，
表示质点i和j的相对运动速度，即lij对时间的导数。上面分析中，i、j表示质点集合M中任意两个相邻的、并由弹簧相连接的质点。
质点i的动态加速度ai可以表示为这里mi是质点i的质量。这里的弹簧静止长度rij是一个随时间和空间分布而变化的变量，模拟叶脉的弹簧由于叶脉中的维管束中拥有较多的纤维，相对更硬一些，因此这些弹簧的静止长度较短，弹簧硬度相对较硬；模拟叶肉细胞的弹簧则相对较软，静止长度一开始相对较长。
4.3下层模型 下层模型表示叶子叶脉的树状结构，而叶脉的结构也正是决定整个模拟过程中最重要的环节。以叶脉节点作为表示叶脉细胞的下层质点，同一根叶脉上每两个相邻的下层质点与其正上方的两个上层质点构成一个矩形，这四个质点彼此之间通过弹簧两两相连，这样就给出了一种矩形中带十字交叉结构的下层质点-弹簧模型，如图9所示，每个交叉结构都有6根弹簧构成，下层模型就主要由这些相互连接的交叉结构组成。在图9中，v1、v2和v3是表示叶肉细胞的上层质点；b1、b2和b3是表示叶脉细胞的下层质点；弹簧sb1v2，sb2v1，sb2v3和sb3v2这4条相交组成“X”状的交叉弹簧，以及sb1b2和sb2b3这两条控制下层叶脉沿着叶脉方向链接的弹簧，用来表示下层较硬的弹簧结构表示叶脉，并控制下层叶脉的结构；同时，弹簧sb1v1，sb2v2和sb3v3这三条竖直方向上上下连接叶脉和叶肉的弹簧是连接上下层结构的，用来控制叶肉和叶脉之间的相对位置和方向关系，保证同一位置的叶脉和叶肉结构之间可以严格地上下邻接。
由于模型中所有弹簧都遵守胡克定律，因此可用方程(3)来进行变形模拟。由于现实世界中，叶脉在不同季节随着水分含量的变化，会导致其物理特性的变化，如硬度、延展力、体积等的变化。因此在下层叶脉模型中，理论上讲，较硬弹簧的弹性系数和阻尼系数是随着时间变化的，并且空间分布也不同，也就是说是随时间-空间变化的，这与上层模型的情况相同。但是，由于水分流失对于叶脉的硬度影响不像对于叶肉细胞的影响那么明显，因此，为了简化整个系统中的参数设置，假定表示叶脉的弹簧和连接上下层模型的弹簧的变形系数都不会发生改变，于是下层模型中这些弹簧的变形系数(包括弹性系数hij、阻尼系数dij和弹簧的静止长度rij)被简化设置为不随时间发生改变。通过本发明进行的一系列测试中，这种简化可以很好地接近真实结果。
图9的右上部分展现了，当上层叶肉结构收缩时，整个叶脉结构是如何向上弯曲的。
图9下半部分展示了沿叶脉纵切面向上卷曲的形状模拟。
4.4上下层模型的结合 在具体仿真过程中，采用质点-弹簧系统实现这种双层模型的模拟，其中代表叶肉细胞的上层模型由一层较软的弹簧相互连接，组成网络，这个代表叶肉细胞的弹簧网就是由之前的三角化得到的叶肉表面的三角网格。而代表叶脉细胞的下层模型由包含交叉结构的较硬弹簧所组成。如图9所示，在每个叶脉片段的位置，将叶脉上方的叶肉抽象成为一条比较柔软的弹簧线(这条线是叶肉弹簧网络中的一个部分)，而下层的叶脉抽象成为一条相对更加刚性的弹簧线，中间由一组交叉结构弹簧将上层的叶肉部分质点与下层的叶脉部分质点进行连接。
这样，两层结构再由弹簧垂直相连，最终构成叶子的整体力学结构，如图10所示，上层结构是一个平面的质点-弹簧网络，下层模型是由每个叶脉处，代表叶脉的较硬的弹簧片段组成的质点-弹簧结构。
5.物理模拟运算 有了上述的能够真实模拟叶子物理结构的质点-弹簧层次结构，就可以开始对叶子的几何造型进行物理模拟，即在本发明的质点弹簧模型结构构建好后，通过Visual C++编写的程序按照时间变化，进行逐步的模拟，随着弹簧的收缩，每经过一个时间段，叶子的构型开始在结构的控制下开始发生变化，并逐步的模拟出自然状态下植物叶子的形态。
具体的物理模拟程序中，本发明编写的程序使用了欧拉迭代的方法对整个质点-弹簧模型的形态变化和运动进行模拟，具体方法将在下面进行阐述，原理可以参见文献[24]。初始状态拥有由前面方法各个步骤生成的质点集合M和弹簧集合S这个质点和弹簧构成的模型，以及它们相关的参数。一般，将最末端叶脉段的宽度设为1，这也就表示最末端叶脉弹簧的硬度是标准1，在叶子质点-弹簧双层模型中也将表示叶肉的上层弹簧硬度设置为标准1，这里的生物学含义就是最末端的叶脉硬度和弹性收缩系数跟叶肉细胞差不多，但是随着叶脉的加粗和维管束和纤维比例的增加，那些高级的、更粗的叶脉的硬度增加、体积变化会更小。如前所述，根据公式推算出整个叶脉结构中各级叶脉段的宽度，由此宽度来定义更上一级别叶脉、乃至主干叶脉弹簧的硬度，一般其硬度在1～50之间。物理模拟过程中具体表现弹簧硬度的参数就是hij，因此，按照以上规则，通过对包含叶肉、叶脉和上下层连接结构的整个质点-弹簧模型的弹性系数进行参数设置。同时，为了简化整个系统中的参数设置，设定表示叶脉的弹簧和连接上下层模型的弹簧的变形系数不随时间发生改变，表示叶肉的弹簧的弹性系数hij和阻尼系数dij也不随时间发生改变，只有表示叶肉的弹簧的静止长度rij(t)随时间变化。
具体的物理模拟步骤如下，假定模拟过程从时间0开始，模拟到了时间t 1.根据由弹簧sij相连接的质点i与j的位置信息，计算出弹簧sij的受力情况fij(t)，这个过程要对弹簧集合S中所有弹簧进行计算。
2.将由sij计算出的弹簧力施加到与其相连的两个质点i与j上，即对于任何质点i，这步将所有弹簧的受力计算出来，施加到质点上。
3.对于质点集合中的每一个质点，例如质点i，由其受力、质量，根据牛顿力学公式计算出其在这个时刻加速度ai(t)。
4.对于质点集合中的每一个质点，例如质点i，其初始的运动速度vi(0)＝0，在t时刻由其加速度ai(t)和速度vi(t)，根据运动学公式，计算出下一个时刻该质点应该处在的位置Pi(t+1)及其新的速度vi(t+1)。
5.更新所有质点的速度、加速度和位置信息到下一个计算时刻(t+1)，并更新表示叶肉的弹簧的变形系数，即通过减小叶肉弹簧的静止长度rij(t)来体现叶肉的收缩。
6.重复步骤1。
不断的重复上面的过程，每个质点的位置和速度也就随着这个物理模拟过程慢慢的发生变化，模拟过程中，模型几何形状的变化，也就是被模拟的叶子几何形态的变化，最终可得到模拟的最终结果，从而也就是叶子老化的最终形态。
本发明通过对于植物细胞生物学和力学结构的分析发现，在叶子枯萎过程中，叶肉和叶脉由于各自不同的组织结构，导致了收缩比例不同，从而使叶子的形态发生变化。在物理模拟中，首先将质点-弹簧模型引入植物模拟领域，用来模拟植物在特定情况和环境下的形态变化。通过建立关于叶子基本结构的双层质点-弹簧模型，并对上下两层的不同参数的合理设置，很好地表现了叶肉和叶脉的不同力学特性。双层模型的相互作用，决定了叶子最终变形的方式和效果，由此可以非常真实地模拟不同种类植物叶子在干枯、老化过程中几何形态的变化。
本发明所提供的植物叶子模拟方法使叶子的绘制结果更加贴近自然真实。由于植物组织拥有相似的生物学构成和物力特性，本方法同时可以扩展到植物整体造型的模拟，对植物的枝干进行基于物理的模拟建模。通过本发明的方法可以在虚拟现实的相关应用中构建更加逼真的植物场景，也可以让植物的真实感绘制在游戏、影视、生物、环保等领域得到进一步应用。
图1是对植物叶子的扫描图像设置不同的分界线取值后获得的不同的二值图像。
图2是叶子轮廓的步进方格法提取方式示意图。
图3是步进方格法提取二值图像时方格四角颜色可能出现的16种情况示意图。
图4是步进方格法中，不同大小的方格所提取的叶子轮廓的细节示意图。
图5是Runions的基于营养点方法生成叶脉的步骤示意图。
图6是现有方法和本发明方法生成的叶脉与真实叶子图像的比较效果图。
图7是拥有叶脉和叶子轮廓后对叶子进行三角网格化后的结果示意图。
图8是细胞在微观情况下力学模型和受力分析示意图，其中A是对单个细胞受力分析的示意图，B是多个细胞受力分析和弹簧模拟情况示意图，C是抽象出来的质点弹簧网络模拟细胞物理结构的示意图。
图9是用于构建叶子下层模型的交叉结构的示意图。
图10是用本发明的双层质点-弹簧模型构建的叶子整体力学结构示意图。
图11是用本发明的方法模拟的爬山虎叶子的变形的绘制结果。
图12是用本发明的方法模拟的椿树叶子和梧桐叶子的变形结果。
图13是用本发明的方法模拟变形的叶子通过光子跟踪方法绘制的结果。
图14是用本发明的方法模拟的多片椿树叶子在树林中绘制的结果。
具体实施例方式 下面通过实施例，结合附图进一步说明本发明，但不以任何方式限制本发明的范围。
根据本发明的基于生物和物理的仿真系统，构造了不同种类的叶子在不同季节中的各种形态。图11是对爬山虎叶子在不同外界条件和参数下，对其构型物理模拟的结果。图12展示了其中的三个实验场景由上到下依次是爬山虎叶子的变形序列、椿树叶子变形的序列和梧桐叶子变形的序列。图13和图14是将变形的叶子模型结果进行全局光照绘制的结果。具体实施过程中同时将不同季节的叶子纹理应用到变形后的模型上。其中，椿树和爬山虎叶子都展示了正面和背面观察的变形效果。由于各种叶子的形状、大小和纹理的不同，绘制结果的成像分辨率定为800×600。所有实验结果均在1.6GHz Athlon CPU，1GB内存和Nvidia 7300显卡的普通PC上，并采用Ray-tracing算法绘制。
大部分叶子结构的质点弹簧模型由1500到5000个质点构成，由于质点弹簧模型在模拟的过程中，通过阻尼振动，逐渐到最后的稳定状态，因此其计算过程通常需要较长时间，具体的模拟时间如下表
参考文献Honda，H.1971.Description of the form of trees by the parameters of the tree-like body.Journal of Theoretical Biology 31，331-338.Weber，J.，and Penn，J.1995.Creation and rendering of realistic trees.In SIGGRAPH’95Proceedings of the 22nd annual conference on Computer graphics and interactivetechniques，ACM Press，New York，NY，USA，119-128.Reeves，W.T.，and Blau，R.1985.Approximate and probabilistic algorithms for shadingand rendering structured particle systems.In SIGGRAPH，313-322.Bloomenthal，J.1985.Modeling the mighty maple.In SIGGRAPH，305-311.Holton，M.1994.Strands，gravity，and botanical tree imagery.Comput.Graph.Forum 13，1，57-67.Fournier，A.，and Grindal，D.A.1986.The stochastic modeling of trees.In Proceedings ofGraphics Interface’86/Vision Interface’86，Canadian Information Processing Society，Toronto，Ont.，Canada，Canada，164-172.Viennot，X.G.，Eyrolles，G.，Janey，N.，and Arqu′es，D.1989.Combinatorial analysis oframified pattems and computer imagery of trees.In SIGGRAPH，31-40.Lintermann，B.，and Deussen，O.1999.Interactive modeling of plants.IEEE Comput.Graph.Appl.19，1，56-65.de Reffye，P.，Edelin，C.，Franc，on，J.，Jaeger，M.，and Puech，C.1988.Plant modelsfaithful to botanical structare and developmentr.In SIGGRAPH，151-158.Prusinkiewicz，P.，and Lindenmayer，A.1990.The algorithmic beauty of plants.Springer-Verlag New York，Inc.，New York，NY，USA.Lefebvre，S.，and Neyret，F.2002.Synthesizing bark.In Rendering Techniques(Eurographics Workshop on Rendering).Federl，P.，and Prusinkiewicz，P.2004.Finite element model of fracture formation ongrowing surfaces.In International Conference on Computational Science，138-145.Dengler，N.，and Kang，J.2001.Vascular patterning and leaf shape.Current opinion inplant biology 4，1，50-56.Judd，W.W.，Campbell，C.S.，Kellogg，E.A.，Stevens，P.F.，and Donoghue，M.J.1999.Plant SystematicsA Phylogenetic Approach.Sinauer Associates，Sunderland，MA.Runions，A.，Fuhrer，M.，Lane，B.，Federl，P.，Rolland-Lagan，A.-G.，and Prusinkiewicz，P.2005.Modeling and visualization of leaf venation pattems.ACM Trans.Graph.24，3，702-711.Sung Min Hong*，Bruce Simpson，Gladimir V.G.Baranoski.Interactive venation-basedleaf shape modeling.In Proceedings of CASA 2005，415-427Murray，C.D.1926.The physiological principle of minimum work.In Proceedings of theNational Academy of Sciences 12，207-214.S.S.Lafon.Diusion Maps and Geometric Harmonics.PhD thesis，Yale University，2004.William E.Lorensen and Harvey E.Cline.Marching cubes A high resolution 3d surfaceconstruction algorithm.In Proceedings of the 14th annual conference on Computergraphics and interactive techniques，pages 163-169.ACM Press，1987.Geoff Wyvill，Craig McPheeters，and Brian Wyvill.Data Structure for Soft Objects.TheVisual Computer，2(4)227-234，February 1986.Brenda S.Baker，Eric Grosse，and C.S.Rafferty.Nonobtuse Triangulation of Polygons.Discrete and Computational Geometry 3(2)147-168，1988.Marshall Bern and David Eppstein. Mesh Generation and Optimal Triangulation.Computing in Euclidean Geometry(Ding-Zhu Du and Frank Hwang，editors)，LectureNotes Series on Computing，volume1，pages 23-90.World Scientific，Singapore，1992.L.Paul Chew.Constrained Delaunay Triangulations.Algorithmica 4(1)97-108，1989.Andrew Witkin and David Baraff.2001，Physically Based Modeling.Course Notes onACM Siggraph 2001 Conference.
权利要求
1.一种植物叶子的模拟方法，包括下列步骤
1)植物叶子轮廓提取对植物叶子进行图像扫描，得到叶子的平面图像，然后对该图像进行二值化，获得二值图像，最后用步进方格法提取叶子轮廓；
2)植物叶脉构造利用基于营养点的叶脉生成方法生成叶脉，所生成的叶脉由连接各叶脉节点的线段组成，在生成叶脉的过程中设定主脉和中脉的生长速度是次级叶脉的1.5-4倍；同时，根据Murray定律对每段叶脉的直径进行设置其中dparent、dchild1和dchild2分别表示叶脉分枝处分枝前的叶脉段的直径和分枝后形成的两个叶脉段的直径，n＝1或2；
3)叶子表面模型三角化以叶脉和叶子轮廓为约束，对叶子表面进行受限制的Delaunay三角化，构建出叶子表面的三角网格结构；
4)构造叶子质点-弹簧双层模型根据步骤3)构建的三角网格结构，以三角网格的顶点作为表示叶肉细胞的质点，三角网格中的边表示连接相邻质点的弹簧，构建出表示叶肉组织的上层质点-弹簧网络模型；根据步骤2)构造的叶脉，以叶脉节点作为表示叶脉细胞的下层质点，同一根叶脉上每两个相邻的下层质点与其正上方的两个上层质点构成一个矩形，这四个质点彼此之间通过弹簧两两相连，构建出矩形中带十字交叉结构的下层质点-弹簧模型，并与上层质点-弹簧网络模型连接起来，从而形成一个可表达叶子力学结构的质点-弹簧双层模型；在该双层模型中，用i、j表示质点集合M中任意两个相邻的、并由弹簧相连接的质点，用sij表示连接质点i和质点j的弹簧，pi和pj分别表示质点i和质点j的位置，Si表示所有与质点i相连的弹簧的集合，则作用在质点i上的力是
其中，lij表示从pi到位置pi的向量；hij(t)表示弹簧sij的弹性系数；dij表示弹簧sij的阻尼系数；rij(t)表示弹簧sij的静止长度；
表示质点i和j的相对运动速度，即lij对时间的导数；
5)进行物理模拟运算根据步骤4)构造的叶子质点-弹簧双层模型编写程序，按照时间变化逐步模拟出自然状态下植物叶子的形态。
2.如权利要求1所述的植物叶子模拟方法，其特征在于所述步骤1)中步进方格法提取叶子轮廓的具体方法是假设叶子的轮廓是封闭和连续的，用方格对叶子二值图像从左上到右下的顺序进行扫描线方式的水平扫描，直到方格的某个角进入叶子内部的区域，然后方格停止扫描线方式的水平扫描，转而开始沿着叶子的边界以逆时针方向运动，直到得到整个叶子的多边形轮廓。
3.如权利要求1所述的植物叶子模拟方法，其特征在于所述步骤2)中设定主脉和中脉的生长速度是次级叶脉的2倍。
4.如权利要求1所述的植物叶子模拟方法，其特征在于所述步骤3)使用CGAL库提供的Exact_predicates_inexact_constructions_kernel类、Delaunay_mesher_2类、Delaunay_mesh_face_base_2类、Constrained_Delaunay_triangulation_2类这4种数据结构来组织以叶子轮廓和叶脉作为限制的线段和点数据，并使用refine_Delaunay_mesh_2函数进行受限制的Delaunay三角化，用数据结构Criteria进行三角形最小角度和最小面积控制。
5.如权利要求4所述的植物叶子模拟方法，其特征在于所述步骤3)中选取参数Criteria(0.005，0.5)，表示三角形最小边的边长不能小于0.005，最小角的角度不能小于0.5弧度。
6.如权利要求1所述的植物叶子模拟方法，其特征在于在步骤5)进行物理模拟运算时设定表示叶脉的弹簧和连接上下层模型的弹簧的弹性系数hij、阻尼系数dij和静止长度rij不随时间发生改变；表示叶肉的弹簧的弹性系数hij和阻尼系数dij也不随时间发生改变。
7.如权利要求6所述的植物叶子模拟方法，其特征在于在步骤2)中将最末端叶脉段的直径设为1，根据公式计算得到其它各级叶脉段的直径，在步骤5)中以叶脉段的直径表征叶脉弹簧的硬度，最末端叶脉弹簧的硬度即为标准1，同时将表示叶肉组织的上层质点-弹簧网络模型中弹簧的硬度也设置为标准1，而弹簧的硬度具体体现为弹性系数hij的大小，由此对整个质点-弹簧双层模型中的弹簧的弹性系数进行参数设置。
8.如权利要求6或7所述的植物叶子模拟方法，其特征在于所述步骤5)进行物理模拟的具体过程是，设定叶子质点-弹簧双层模型的初始参数后，从t＝0开始依次重复下列的步骤a～e
a.对于每个弹簧sij，根据由弹簧sij相连接的质点i与j的位置信息，计算出在时刻t弹簧sij的受力情况fij(t)；
b.将由sij计算出的弹簧力施加到与其相连的两个质点i与j上，即对于任何质点i，
c.对于每一个质点i，由其在时刻t的受力fi(t)和质量mi，根据牛顿力学公式计算出其在这个时刻的加速度ai(t)；
d.对于每一个质点i，其初始的运动速度vi(0)＝0，在t时刻由其加速度ai(t)和速度vi(t)，根据运动学公式，计算出下一个时刻t+1时该质点应该处在的位置pi(t+1)及其新的速度vi(t+1)；
e.更新所有质点的速度、加速度和位置信息到下一个计算时刻(t+1)，并更新表示叶肉的弹簧的的静止长度rij(t)。
全文摘要
本发明提供了一种基于物理的对植物叶子进行模拟的方法，主要包括植物叶子轮廓提取、植物叶脉构造、叶子表面模型三角化、构造叶子质点—弹簧双层模型结构和进行物理模拟运算5个步骤。本发明基于植物学和物理原理，建立了质点—弹簧双层模型来表达叶子的力学结构，对上下两层的不同参数的合理设置，很好地表现了叶肉和叶脉的不同力学特性，可以非常真实地模拟不同种类植物叶子在干枯、老化过程中几何形态的变化。由于植物组织拥有相似的生物学构成和物理特性，本方法同时可以扩展到植物整体造型的模拟，对植物的枝干进行基于物理的模拟建模。
文档编号G06T17/00GK101251932SQ200810102119
公开日2008年8月27日 申请日期2008年3月18日 优先权日2008年3月18日
发明者吴恩华, 迟小羽, 斌 盛, 王文成 申请人:中国科学院软件研究所