沉积物作为水体环境中重要的组成部分，既是各种污染物的汇集地，也是对水体水质具有潜在影响的次生污染源.沉积物是具有复杂结构的多孔介质(dullien，1992)，其孔隙结构具有非均质、不连续和复杂不规则的结构特点，孔隙体之间通过喉道相互连接，构成的空间网络体系是物质迁移运输过程得以发生的前提条件，污染物在其中的宏观迁移和扩散过程受微观的孔隙度、孔隙连通性和曲率等系统形态及其几何尺寸的影响(brakel et al., 1974)，在相同孔隙度的多孔介质中，曲率增大或孤立孔数量增加都会导致物质扩散系数和渗透系数的减小，因此，要定量地研究沉积物孔隙结构对污染物运动过程的影响机制，必须揭示沉积物的真实三维物理结构并对其进行准确描述.
近几十年来计算机图像处理技术与x射线断层摄影技术(computed tomography，ct)的快速发展，为研究多孔介质孔隙结构提供了有力的工具(elliot et al., 2007)，它是目前最先进的无损检测手段，具有图像清晰直观、密度分辨率高、图像数字化等优点(高丽娜等，2009).ct最早的应用是在医学诊断领域，20世纪70年代中后期ct技术在工业上得到发展，并成功应用于土壤微观结构的研究中，随后迅速发展，在各个学科领域得到广泛的应用.yoshito nakashima根据ct扫描结果，编写mathematica程序分析得到了多孔岩石孔隙在三维空间上的连通性.温胜芳(2013)从沉积物采样制样方法上进行创新，克服了未成岩沉积物含水量高、结构松散、易扰动等问题，并构建了一整套未成岩沉积物三维物理结构研究方法，弥补了该领域研究方法的不足，为研究分析沉积物中微观孔隙结构提供了基础.
巢湖位于安徽省中部，是我国五大淡水湖之一，河流输入性和浅水特点导致巢湖沉积物沉积环境复杂.由于流域内不同区域的水土流失、土地利用等有很大的差异，导致巢湖湖区存在不同的沉积物组成和污染分区(温胜芳，2013)，这些特征为沉积物结构及其环境效应研究提供了很好的差异性.因此，本文选择巢湖西、中、东3个湖区6个采样点的沉积物作为研究对象，利用ct技术，计算分析各采样点表层沉积物孔隙结构参数特征及垂直变化特点，揭示巢湖沉积物真实三维孔隙物理结构，以弥补模型构建方法对于认识沉积物真实孔隙结构研究的不足.
2 材料与方法
2.1 研究区域概况
巢湖(31°25′~31°43′ n，117°16′~117°50′ e)位于安徽省中部，江淮丘陵之间，湖岸线总长184.7 km，东西长61.7 km，南北宽20.8 km，最大水深3.77 m，平均水深2.69 m，现有水域面积大约为770 km2，蓄水量20.7×108 m3(王苏民等，1998).流域总面积达16659 km2，以巢湖与裕溪河之间的闸门分为两个部分(余秀娟，2012)，属亚热带与暖温带过渡性的副热带季风气候区，气候温和湿润，年平均温度在15~16 ℃之间，多年平均降水量为1000~1158 mm.流域内地形地貌西高东低、中间低洼平坦(顾成军等，2005).流域内共有河流33条，分别属于南淝河-店埠河、杭埠河-丰乐河、白石山河、柘皋河、派河、裕溪河等7条水系，从南、西、北呈放射状直接或间接汇入巢湖，再由裕溪河注入长江(高超等，2009).1962年建成巢湖闸控制水位后，基本上切断了巢湖与长江在自然状态下的水量交换，巢湖由自然状态下的过流性湖泊，成为半吞吐型过水湖泊.由于入湖上游坡陡流急，挟带大量泥沙汇入湖内沉积，多年来，巢湖因流域河流输沙影响而导致淤积严重(阎伍玖等，1998).
2.2 样品点及采样方法
根据巢湖泥沙沉积及污染情况，分别在西、中、东3个湖区内各设置采样点2个(图 1)，6个采样点分别对应西部污染区(w1)、西部泥沙较重区(w2)、中部泥沙较轻区(m3)、中部泥沙区(m4)、东部泥沙较轻区(e5)和东部水质较好区(e6).2014年7月，依托自重力沉积柱采样平台，自制沉积物柱采样管，利用液氮低沸点大热容的特点，解决常规采样方法对流态沉积物结构的扰动问题，原位冷冻固定沉积物，获得6根原位沉积物柱状样，置于装有液氮的保温壶中保存，其中，m3点样品由于保存问题，无法进行后续处理，导致相应数据丢失.
图1 巢湖湖区采样布点示意图
2.3 孔隙结构研究方法
借鉴温胜芳(2013)构建的整套包含样品制备方法、数据获取方法、图像解译方法、参数计算方法的未成岩沉积物三维物理结构研究方法，计算分析巢湖沉积物孔隙结构参数，具体过程如下.
样品制备：通过初步切割-冷冻干燥-lr white树脂包埋-切割-增加参比薄片(树脂和矿物质)等步骤，将冷冻沉积柱样品(φ2.6 cm×10 cm)制备成 5 mm×5 mm×8 mm 的稳定的包括树脂和矿物质参比的满足ct扫描和后续图像解译要求的小体积样品.实验中通过选用低粘度的lr white树脂及包埋过程中各操作的控制，保证样品原始孔隙结构.
数据获取：选择工业锥束显微ct系统(中国科学院高能物理所90 kev cone beam micro-ct)，在70 kev的电压、90 μa的电流条件下对样品进行扫描，利用fdk滤波反投影算法，重建被测物体的三维图像，并对几何偏差进行校正，得到样品16位无符号整型图像序列.ct图像空间分辨率达到6 μm.
图像解译方法：根据实际样品和标准样品灰度分布图像序列，利用vg图像处理软件确定分割阈值，并将16位无字符数据转换为8位无字符数据，输出兴趣区的tiff图像序列用于后续分析.
参数计算方法：用mathematica程序对三维孔隙进行连通性标记，并在标记的连通孔内用粒子随机行走模型计算各向异性曲率等参数，共包含3个程序：trim程序用于对图像进行兴趣范围选择;clabel程序用于对孔隙群进行标记，通过输入原始灰度图像序列、孔隙分割阈值和孔隙群等级颜色文件，最终输出一个包含各孔隙群的体积、表面积、中心点坐标、是否临近边界等参数的文件，分析孔隙连通性和孔隙度等特性;rwalk程序用于在标记的孔隙群的连通孔隙中进行粒子随机行走模拟，设定随机行走的粒子数量和行走时间，最终根据均方位移与时间的斜率进行曲率的估算.
3 结果与分析
3.1 不同湖区表层沉积物孔隙结构特征
对不同湖区表层沉积物孔隙参数进行计算，各点计算结果如表 1所示.东部湖区e6点连通状况最好，孔隙度为0.75，连通孔体积达2.53 mm3，孤立孔数量最少，中部m4点连通性较差，孔隙度仅为0.28，其他点的孔隙度变化范围为0.51~0.70，连通孔体积在1.71~2.37 mm3之间变化，其中，西部湖区w1点无论是像素大于4(864 μm3)的相对较大孤立孔，还是像素小于4的孤立孔，均明显高于其他点，表明了巢湖表层沉积物孔隙连通性较高但细小孔隙多的特点，有利于沉积物中物质的转化和跨界面迁移;连通孔s/v(表面积与体积之比)在不同采样点差别较大，变化范围为159~296.
表1 巢湖不同湖区表层沉积物孔隙连通性(沉积物体积3.375 mm3，250×250×250像素)
通过mathematics软件显示沉积物各点三维孔隙结构，如图 2所示，图中深蓝色部分为沉积物基质，黄色部分为连通孔，浅蓝部分为孤立孔.
图2 巢湖不同湖区表层沉积物的孔隙结构标记图
根据以上结果，利用rwalk程序对粒子在自由空间与沉积物孔隙中的扩散进行随机行走模拟，结果见图 3(以w1点1~2 cm深度内样品为例)，其中，空心方框和三角分别表示行走轨迹的起点和终点.由于沉积物结构的各向异性，根据不同方向的均方位移与时间统计结果可以得到不同方向的曲率值.
图3 巢湖w1样点表层沉积物连通孔中粒子随机行走结果
巢湖不同湖区样点表层沉积物(1~2 cm深度内)的曲率见表 2，结果表明，巢湖沉积物曲率较大，且各向具有差异性.m4点综合曲率最高，为9.31，说明沉积物对离子扩散阻碍作用最强，其他几个点的曲率在2.83~6.79之间变化.此外，不同方向上的曲率值不同，w1、w2点在z方向的曲率较大，m4和e5、e6点分别以x和y方向的曲率相对较大，反映了不同沉积物结构在空间各向上的不同.
表2 巢湖不同湖区表层沉积物曲率
3.2 沉积物孔隙结构垂直变化特征
对巢湖西部湖区w1点不同深度沉积物孔隙参数分布特征进行分析，孔隙三维结构见图 4，参数计算结果如表 3所示.孔隙连通性的剖面数据显示，w1点沉积物剖面的孔隙度与连通孔体积随深度增加而减小，0~10 cm内变化不大，分别由0.51、1.71 mm3降低到0.44、1.48 mm3，10 cm后变化幅度明显增大，17 cm后孔隙度为0.01，说明沉积物与上覆水物质的交换主要集中在沉积物上层，但深度变化对s/v影响较小，变化范围为272~299，变化幅度明显低于其表层沉积物在不同湖区的变化;孤立孔数量、体积与表面积随深度增加均呈现出减小的趋势，但在13~14 cm深度内出现异常，表明沉积物孔隙度等参数除受到沉积深度影响外，还可能受到当时底栖动物活动及沉积条件的共同作用.
图4 w1点不同深度沉积物的孔隙群标记图
表3 巢湖w1点不同深度沉积物孔隙群统计表(3.375 mm3，250×250×250像素)
通过mathematics软件显示沉积物不同深度样品三维孔隙结构，如图 4所示，图中深蓝色部分为沉积物基质，黄色部分为连通孔，浅蓝部分为孤立孔.随后利用rwalk程序对粒子在自由空间与沉积物孔隙中的扩散进行随机行走模拟，计算曲率.
巢湖w1点不同深度样品的曲率见表 4，其中，17~18 cm深度范围内孔隙度极低，仅为0.01，因此，未进行随机模拟行走计算曲率.从表 4可以看出，随深度的变化，综合曲率变化范围为5.89~8.03，且均表现为z方向上曲率相对较大，表明w1点孔隙水中污染物向上的扩散作用会受到较大的阻碍.
表4 巢湖w1点不同深度沉积物曲率
4 结论
本文利用已有的沉积物原位采样方法及孔隙结构研究方法，通过ct技术及相关软件，得到了巢湖沉积物真实的三维孔隙结构，并计算得到了巢湖不同湖区采样点沉积物的孔隙度、孔隙连通性、曲率等参数，分析了巢湖各点沉积物结构特征及西部湖区沉积物结构参数垂直变化特点，主要结论如下：
1)各湖区表层沉积物孔隙度不同，东部湖区e6点连通状况最好，孔隙度为0.75，连通体积达2.53 mm3;中部m4点连通性较差，孔隙度仅为0.28，曲率在2.83~9.31内变化，且在各个方向上表现出明显的差异性;w1、w2点在z方向的曲率较大，m4、e5、e6点以x和y方向的曲率相对较大，反映了不同沉积物结构在空间各向上的不同，进而会导致沉积物中污染物在空间各个方向上的运输与扩散的差异.此外，孤立孔体积、连通性等参数也有所差别.除m4点外，巢湖表层沉积物较大的孔隙度和孔隙连通性及较低的综合曲率均有利于沉积物中物质与上覆水体的迁移运输.
2)w1点沉积柱在垂直剖的面数据显示，孔隙度随深度增加而减小，0~10 cm内变化不大，由0.51降低到0.44，10 cm后变化幅度明显增大，17~18 cm深度内孔隙度为0.01，说明沉积物与上覆水物质的交换主要集中在沉积物上层，但综合曲率变化不大，说明这种情况下，曲率不是限制下层沉积物中物质迁移运输的主要因素.此外，深度变化对s/v影响较小，变化范围为272~299，变化幅度明显低于其表层沉积物在不同湖区的变化，孤立孔数量、体积与表面积随深度增加均呈现出减小的趋势，但在13~14 cm深度内出现异常，这可能受到当时底栖动物活动及沉积条件的影响.
3)本文通过应用ct技术构建了巢湖沉积物三维孔隙结构，直接计算得到fick 第一定律通量估算中的沉积物孔隙度φ值，并可以减少通过经验公式计算分子扩散系数 ds的误差，从而提高了孔隙水扩散模型定量估算沉积物内源负荷的准确性.