城市表层土壤重金属污染分析大学生数学建模一等奖论文.doc

城市表层土壤重金属污染分析摘要 本文对城市表层土壤中重金属的污染进行分析，得出土壤中重金属污染的主要原因，污染的传播特性以及各污染源的位置，这些分析结果对于治理土壤重金属污染有很大帮助。 对于问题一，把取样点的位置坐标分别作为x轴和y轴，各种元素的浓度作为z轴，对各城区每种金属的浓度分别进行三次多项式插值，得到浓度与取样点位置的函数关系，并画出8种重金属元素的浓度在该城区的空间分布图以及单个元素浓度的空间分布图，如图1-图10。用单因子污染指数法和内梅罗综合污染指数法相结合对该城区内不同区域重金属的污染程度进行分析，依据内梅罗指数土壤污染评价标准对不同区域重金属的污染程度进行分级，最终得出每个区内各种金属元素的污染等级，结果见表1。对于问题二，首先建立了综合权重模型，计算了各区内8种金属污染物的综合权重并对其进行排序，得出金属污染最为严重的区为工业区和交通区，其中工业区的综合权重达到105.76，交通区的综合权重为52.62。其次，运用主成分分析法，提取出每个区的主成分，然后根据每个主成分与8种重金属污染物的相关程度来确定出该区的主要重金属污染物。其中，一区的主要重金属污染物为Pb，Zn，Ni；二区为Cr，As；三区Zn，As，Hg；四区Cu，Hg；五区的主要重金属污染物是Cd，Pb，Hg。对于问题三，根据三个空间方向上8种重金属污染物浓度的聚类点趋势图分析，建立非线性回归模型，用最小二乘法求解得到每种金属浓度的传播特性曲线函数，从而确定出污染源的位置，结果见表6和表7。对于问题四，对问题三中模型的优缺点进行了评价与分析。由于降雨量、风速、人流量、交通流量也是重金属污染的影响因素，因此综合考虑这些因素建立积分方程模型分析重金属污染物对地质环境的影响。 关键词：污染等级 综合权重 主成分分析 曲线拟合 一 问题重述随着城市经济的快速发展和城市人口的不断增加，人类活动对城市环境质量的影响日显突出。对城市土壤地质环境异常的查证，以及如何应用查证获得的海量数据资料开展城市环境质量评价，研究人类活动影响下城市地质环境的演变模式，日益成为人们关注的焦点。城区按照功能分为生活区、工业区、山区、主干道路区及公园绿地区等，分别记为1类区、2类区、5类区，不同的区域环境受人类活动影响的程度不同。根据附件1中列出的采样点的位置、海拔高度及其所属功能区等信息，附件2中列出的8种主要重金属元素在采样点处的浓度，以及附件3中列出的8种主要重金属元素的背景值来对某城市城区土壤地质环境进行调查、研究。并通过数学建模来完成下面任务：(1) 给出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度。(2) 通过数据分析，说明重金属污染的主要原因。(3) 分析重金属污染物的传播特征，由此建立模型，确定污染源的位置。(4) 分析第三问中建立模型的优缺点，分析为更好地研究城市地质环境的演变模式，还应收集什么信息？并针对这些信息，考虑如何建立模型解决问题？ 二 问题分析 对于问题一，要求画出8种主要重金属元素在该城区的空间分布，并分析该城区内不同区域重金属的污染程度，对于空间分布，可以运用matlab中的mesh()函数画出8种主要重金属元素在该城区的空间分布图以及各种元素浓度的空间分布图。在分析该城区内不同区域重金属的污染程度时，可以运用单因子污染指数法和内梅罗综合污染指数法相结合对该城区内不同区域重金属的污染程度进行分析，并根据内梅罗指数土壤污染评价标准对不同区域重金属的污染程度进行分级。对于问题二，要求通过数据分析，说明重金属污染的主要原因。从两个方面来分析，一方面，从每个区对所有金属污染物的影响出发，来研究哪些区对重金属污染物的影响比较大，分别计算出各个区内8种金属污染物的综合权重并对其进行排序，找出污染最严重的区；另一方面，用主成分分析法和因子分析法，首先用SPSS软件进行统计分析，得出每个区的主要重金属污染物，然后再根据各个污染物所在区的环境特征及人类活动情况来分析得出重金属污染的主要原因。对于问题三，通过曲线拟合分别得到在x，y，z方向上8种金属污染物浓度的影响曲线，来分析x与y分量对各金属污染物浓度的传播特性的影响及z分量对各金属污染物浓度的传播特性的影响，通过聚类点趋势图分析，得到空间坐标的三个分量与污染浓度的关系，求解最终得到金属浓度的空间曲线函数来确定出污染源的位置。对于问题四，首先对第三问所建立模型进行评价，然后为更好地研究城市地质环境的演变模式，可以收集该城区降雨量、风速、人流量、交通流量等信息，根据这些信息来建立多因素与地质环境演变趋势的模型，进一步研究城市地质环境的演变模式。三 问题假设（1）假设该城区的各项环境指标比较稳定，不会发生突变；（2）假设该城市在一定时间内不会发生重大的自然灾害；（3）假设题目中的所给的各个数据都能够真实的反应实际情况；四 符号说明符 号符 号 解 释功能区第种金属元素的浓度第种金属元素背景值范围的最大值功能区第种金属元素实际测量的总个数横坐标的值纵坐标的值海拔时间段重金属浓度该城区某时刻某地域的风速该城区某时间段某地域的降雨量该城区某时间段某地域的人流量该城区某时间段某地域的交通流量该城区某地域的人流密集度该城区某地域的交通密集度五 模型建立与求解5.1 问题一5.1.1 模型的建立与求解首先要求画出8种主要重金属元素在该城区的空间分布图，主要运用MATLAB绘图来解决，通过使用xlsread()函数将Excel表格中的数据导入到MATLAB中，把取样点的位置坐标分别作为x轴，y轴，各种元素的浓度为z轴，先利用MATLAB的greddata( )函数插值来使得得出的三维图形更加平滑直观，再利用mesh( )函数画出8种主要重金属元素在该城区的空间分布图以及各种元素浓度的空间分布图。5图1 8种重金属元素在该城区的空间分布 图2 该城区5个功能区的空间分布图3 As元素的浓度空间分布图 图4 Cd元素的浓度空间分布图图5 Cr元素的浓度空间分布图图6 Cu元素的浓度空间分布图图7 Hg元素的浓度空间分布图图8 Ni元素的浓度空间分布图图9 Pb元素的浓度空间分布图 图10 Zn元素的浓度空间分布图由图1-图9，可直观的看出重金属As 元素在交通区、工业区的污染程度最重，重金属Cd 元素在交通区、工业区污染程度最重；Cr 元素在交通区、生活区的污染程度最重；Cu 元素在工业区、交通区的污染程度最重；Hg 元素在工业区、交通区的污染程度最重；Ni 元素在在交通区、山区的污染程度较为严重；Pb 元素在交通区、生活区污染程度较重；Zn 元素在各区的工业区、交通区污染程度较重。5.1.2 模型的建立与求解（一）模型的建立要分析该城区内不同区域重金属的污染程度，用单因子污染指数法和内梅罗综合污染指数法相结合进行评价，内梅罗综合污染指数计算公式为： （1）其中为土壤重金属综合污染指数，为土壤中第j种污染物的实测含量，为土壤中第j种污染物的背景值，为土壤中各重金属的平均单相污染指数，为土壤中各金属的最大单相污染指数，内梅罗指数反映了各污染物对土壤的作用，同时突出了高浓度污染物对土壤环境质量的影响，可按内梅罗污染指数，划定污染等级，对8种重金属元素的污染进行划分等级，进而分析并分析该城区内不同区域重金属的污染程度。内梅罗指数土壤污染评价标准如表2。（二）模型的求解 根据附件2和附件3中数据根用内梅罗综合污染指数计算公式对该城区内不同区域重金属的污染程度进行分级，得出各区不同金属污染指数等级表如表1。表1 各区不同金属污染指数等级表各金属AsCdCrCuHgNiPbZn一区1.714.0310.798.797.731.347.8521.16污染等级轻污染重污染重污染重污染重污染轻污染重污染重污染二区3.024.324.1987.75187.391.647.3112.03污染等级重污染重污染重污染重污染重污染轻污染重污染重污染三区1.531.622.562.472.922.691.971.76污染等级轻污染轻污染中污染中污染中污染中污染轻污染轻污染四区4.026.1813.3147.36221.925.103.1627.47污染等级重污染重污染重污染重污染重污染重污染重污染重污染五区1.743.951.535.0818.641.173.8710.19污染等级轻污染重污染轻污染重污染重污染轻污染重污染重污染由表1可得出结论：一区内金属Cd、Cr、Cu、Hg、Pb、Zn为重度污染，金属Ni、As为轻度污染；二区内金属As、Cd、Cr、Cu、Hg、Pb、Zn为重度污染，金属Ni为轻度污染；三区内金属Ni、Cu、Hg、Cr为中度污染，金属As、Cd、Pb、Zn为轻度污染；四区内金属As、Cd、Cr、Cu、Hg、Ni、Pb、Zn为重度污染；五区内金属Zn、Cu、Hg、Pb、Cd为重度污染，金属As、Cr、Ni为轻度污染。表2 土壤重金属内梅罗污染综合指数分级标准等级内梅罗污染指数污染等级P0.7安全（清洁）0.7P1警戒线（尚清洁）1P2轻度污染2P3中度污染 P3重度污染根据表1中的污然指数等级表画出饼图如图11，由饼图可直观的看出各重金属在每个区的污染程度。15图a 功能区一中金属污染程度图b 功能区一中金属污染程度 图c 功能区一中金属污染程度 图d 功能区一中金属污染程度图e 功能区一中金属污染程度图11 各功能区金属污染程度由图11可以直观的看出每个功能区8种重金属元素的污染程度，功能区一中Zn元素污染程度最大，功能区二中Hg元素污染程度最大，功能区三中Hg元素污染程度最大，功能区四中Hg元素污染程度最大，功能区五中Hg元素污染程度最大。5.2 问题二5.2.1 模型的建立要求说明重金属污染的主要原因，首先考虑分别计算出各个区域内8种金属的综合权重并对其进行排序，从而得出对金属污染影响严重的区域，在根据相应区域的环境特征及人类活动情况等来分析得出重金属污染的主要原因。综合权重计算公式为： （2）其中：表示功能区第种金属元素的平均相对污染程度； ：表示功能区第种金属元素的浓度相对于其背景值最大值的超标率；：表示功能区第种金属元素的相对污染程度；其次，用主成分分析法和因子分析法，首先用SPSS软件进行统计分析，得出每个区的主要重金属污染物，然后再根据各个污染物所在区的环境特征及人类活动情况来分析得出重金属污染的主要原因。5.2.2 模型的分析与求解5.2.2.1 模型一的分析与求解从每个区对所有金属污染物的影响出发，来研究哪些区对重金属污染物的影响比较大，分别计算出各个区内8种金属污染物的综合权重并对其进行排序，从而得出对重金属污染的影响最为严重的区域。根据综合权重的计算公式：对附件中的数据进行处理分析得出表3。表3 各个区域内8种金属的综合权重功能区一AsCdCrCuHgNiPbZn平均值1.453.932.918.316.510.394.9810.64超标率0.680.730.450.770.570.270.590.61综合权重24.87功能区二AsCdCrCuHgNiPbZn平均值2.766.951.4829.8474.650.898.5513.12超标率0.610.860.420.940.690.420.810.78综合权重105.76功能区三AsCdCrCuHgNiPbZn平均值1.222.880.9311.1627.240.673.484.91超标率0.170.260.230.260.210.180.180.14综合权重11.23功能区四AsCdCrCuHgNiPbZn平均值1.126.121.7911.8150.710.423.9010.89超标率0.500.780.340.830.520.190.670.72综合权重52.62功能区五AsCdCrCuHgNiPbZn平均值1.423.870.383.109.360.193.675.11超标率0.690.540.200.690.510.140.400.49综合权重14.07根据表3可以看出二区和四区对重金属污染的综合权重比较大，即工业区和交通区是造成重金属污染的主要因素。由于工业区，存在大量的化工厂，其排出的污水和废气中含有大量的重金属污染物，而交通区车流量大并且种类繁杂，由于使用含Pb等重金属污染物的汽油使得排放的尾气中存在严重的重金属污染，因此在工业区和交通区的重金属污染比较严重。5.2.2.2 模型二的分析与求解对单个重金属污染物进行考虑，运用主成分分析法和因子分析法，首先对一区的数据进行分析，步骤如下：（1）将一区的数据导入SPSS中，利用降维的思想，指定八种重金属元素为参与分析的变量，把8个指标转换成较少的几个不相关的综合指标，来进行分析五个区内的主要金属污染物。（2）通过对SPSS软件使用默认值进行因子分析，得出各成分的公因子方差表和主成分分析的因子载荷阵。对得出的表格进行处理的表4和表5。表4 各成分的公因子方差表成分序号各成分特征值各因子特征值占特征值总和的百分比 %各因子方差占总方差的百分比的累积百分比 %13.61645.19945.19921.13314.16559.36531.07513.43272.79740.80710.08382.8850.5246.54489.42460.4455.55694.9870.2362.95497.93480.1652.066100表5 主成分分析的因子载荷比成分12345678As(ug/g)0.699-0.646-0.010.130.17-0.1730.2460.164Cd(ng/g)0.7840.171-0.417-0.23-0.069-0.134-0.2550.204Cu(ug/g)0.6430.2340.493-0.091-0.4350.2810.0750.08Cr(ug/g)0.729-0.2460.024-0.2680.410.405-0.041-0.045Hg(ng/g)0.4920.13-0.4370.704-0.0540.223-0.009-0.036Ni(ug/g)0.686-0.2530.5230.251-0.002-0.22-0.247-0.137Pb(ug/g)0.8030.112-0.348-0.289-0.162-0.1560.15-0.255Zn(ug/g)0.5010.6910.2670.120.363-0.1770.140.058由于确定选取几个成分作为主成分的判定方法有两种：a. 取所有特征值大于1的成分作为主成分；b. 根据累计贡献率达到的百分比值确定（取累计贡献率达到80%）。根据表一首先可以得出特征值大于一的成分序号为1、2、3，其次由表可看出前四个因子的特征值之和占总方差的82.88%，综合两种判别方法最终确定提取三个主成分，使用这些成分相当大的程度上减少了原始数据的复杂性，再根据表二中每个主成分与八种重金属元素因子的相关程度大小进行选择，得出第一主成分与Pb的相关度较高，第二主成分主要与Zn的相关度较高，第三主成分主要与Zn的相关度较高，最终确定出一区的主要重金属污染物是Pb，Zn，Ni。根据上述方法，用SPSS软件我们得出二、三、四、五区的相关数据见附表一至附表八。根据附表一至附表八中数据，用同样的判别方法得出：二区的主要重金属污染物是Cr，As，三区的主要重金属污染物是Zn，As，Hg，四区的主要重金属污染物是Cu，Hg，五区的主要重金属污染物是：Cd，Pb，Hg。总结上述结果可以得到下面结论：对于一区，由于人口密度大，机动车辆多，因此各种垃圾的堆积以及车辆尾气的排放使得此区Pb、Zn、Ni的含量高于其背景值，污染比较严重；对于二区，由于存在大量化工厂，排放的污水及废气中含有大量的重金属污染物特别是Cr，As，因此使得在该区内这些重金属污染比较严重；对于三区，其较其他区而言分布着大量的田地，由于农药、化肥、塑料薄膜的使用，导致土壤中的Zn，As，Hg等含量较多；对于四区由于车流量大，汽车尾气的排放以及汽车轮胎磨损使得Cu，Hg的在交通区的污染比较严重；对于五区，由于人流量很大，且邻近交通区交通发达，仍存在含重金属汽车尾气的排放以及人群的流动使得Cd，Pb，Hg在公园区的含量也高于其背景值。5.3 问题三5.3.1 模型的建立 （一） 模型假设（1）金属污染的浓度只考虑与位置坐标的关系，不考虑气候、水流等其他因素；（2）金属污染源处于较稳定平衡状态，它的位置由其浓度与浓度变化率的两种关系来确定，即、的值。 （二）数据处理将附表二所给各种元素的浓度值减去其背景值的最大值得出该元素是否对环境造成污染，并将为负值的点剔除，即无污染点。（三）模型建立根据处理后的数据我们对重金属污染物的传播特征进行分析，并由此建立模型，来确定污染源的位置。由于空间上金属浓度的变化没确定的规律，因此分别对8个金属污染物浓度在x，y，z方向上的变化曲线进行模拟。由拟合出的曲线可分析出x与y分量对各金属污染物浓度的传播特性是随着其各自的增大近似呈成正弦波波动，z分量对各金属污染物浓度的传播特性是随着其增大呈递减趋势。通过观察拟合图像对数据大致走向进行分析，得到空间坐标的三个分量与污染浓度的关系，即为，通过对上述三个方程进行积分求解，最终得到各金属浓度的空间曲线的函数表达式为： （3）其中：函数各项的系数；：第j种金属浓度的传播特性曲线函数；：第j种金属浓度的传播特性曲线函数对的偏导数；：第j种金属浓度的传播特性曲线函数对的偏导数；：第j种金属浓度的传播特性曲线函数对的偏导数；：表示第j种金属函数在 空间某点的梯度大小； 第j种金属函数在空间某点的浓度与其梯度比值的大小； 第j种金属函数在空间某点的浓度与其梯度乘积的大小；5.3.2 模型的分析与求解根据空间坐标的三个分量与污染浓度的关系 分别对8个金属污染物浓度在x，y，z方向上的变化曲线进行模拟。图11 在x方向上As金属污染浓度的传播特性曲线图12 在y方向上As金属污染浓度的传播特性曲线图13 在z方向上As金属污染浓度的传播特性曲线通过最小二乘法进行曲线拟合，As的传播特性函数为： （4）其曲线图如图14中的红色曲线。 图14 As模拟曲线与对应原始曲线对比图用同样的方法得出其他七种元素的传播特征函数及其相应的传播特性曲线。 对于污染源的位置确定，考虑到如果该污染源为新增污染源，它的浓度变化率应较大，如果该污染源已持续了一段时间则其变化率比较慢，根据这种情况我们将梯度与浓度共同考虑，通过带入已知样本点计算、找出其最大值，它们对应的坐标即为污染源的位置见表6、表7。表6 各金属污染源的位置及其传播方向金属元素坐标功能区传播方向As(19072,8519,36)52.793（0.274i,0.111j,0.329k)Cd(10142,1662,8)43.881(28.698i,-6.023j,-14.305k)Cr(1647,2728,6)17.801(12.574i,-18.406j,-12.027k)Cu（6869,7286,18）10.562(14.771i,-2.224j,-6.013k)Hg(13765,1353,15)9,704(166.400i,-428.721j,353.583k)Ni(22193,12185,79)24.531(-2.024i,-1.307j,-2.191k)Pb(9090,5365,20)59.704(-0.227i,-7.786j,-10.023k)Zn(13797,9621,18)24.941（12.043i,39.720j,61.111k）表7 各金属污染源的位置及其传播方向金属元素坐标功能区传播方向As（7776 ,10613,9）4.944（0.287i,0.053j,-0.441k)Cd（6534,5641,6）20022.572(29.478i, -26.978j,6.076k)Cr（1049,2127,12）4284.884(16.091i,-6.417j,-12.578k)Cu（1049,2127,12）8782.795(11.863i,26.076j,10.061k)Hg（1321,1791,28）677416.64(-189.831i,-405.718j,323.7k)Ni（11646,9381,14）133.432(2.737i,-0.717j,-1.838k)Pb（1049,2127,12）2023.694(-4.415i,-6.174j,-1.955k)Zn（0,1787,4）37308.042(1.658i,30.092j,42.465k)5.4 问题四5.4.1模型的建立与求解对于问题三中的模型，优点：（1）能够得到每种重金属物在该城区内任意一点沿某一方向传播时其浓度变化快慢；（2）从得出的非线性方程我们可以计算出任何一点的金属污染的浓度值；（3）拟合所得到的曲线与原曲线具有良好的相关性。缺点：（1）在模型求解时，数据处理比较繁琐，且需要编写大量程序，最终得出的图表也比较多； （2）在建立该模型时没有考虑多方面因素，例如没有把风向、水流、时间等因素考虑在内，建立的模型不是很严谨； （3）所得到的函数在浓度较大处拟合的不是很完美，误差较大。由上分析可知问题三中模型所得曲线的拟合程度不完美，说明在建模时考虑的影响因素太少使数值偏差较大，那么应该考虑更多可以便捷得到数据的因素，使模型更接近真实值。5.4.2模型的建立 模型假设：（1）该城区较长时间内不会发生自然灾害使地质环境发生重大变化；（2）测量风速时风的方向不发生改变；（3）该城区在某地域某时间段测量的人流量、交通流量具有较高准确性；（4）影响重金属浓度的因素只考虑地理位置、降雨量、人流密集度、交通密集度、风速；为了更好地研究城市地质环境的演变模式，我们考虑到还可以收集的信息包括风速、气候、水源位置、污水排放、工业垃圾、生活垃圾、农药喷洒、人流密集度、交通密集度等因素，各因素与金属浓度的联关系图如下图；易测量、统计或计算因素金属浓度气候（降雨量）水源位置污水排放工业垃圾生活垃圾农药喷洒人流密集度交通密集度风速图16 各因素与金属浓度的联关系图其中可以方便收集到降雨量、风速、人流量、交通流量的数据信息，有了这些信息通过分析我们可以建立模型求解问题。可以使用仪器测量该城区某地域不同时刻的风速，并对该区域中时间段内人流量、交通流量进行统计，再计算出人流密集度、交通密集度，同时测量出该时间段内的降雨量，然后在该地域采样得到重金属的浓度，通过对、模拟曲线的分析分别得到能反映在、分量上传播特性的函数、，因为这些因素之间互不相关，则可计算得到最终的金属浓度的函数表达式： （5）其中，六 模型的优缺点61 模型优点（1）用MATLAB编程时，程序比较简单，并且能够直接调用表格中的数据，画图时可以将五个功能区用不同的颜色在空间图中表现出来，比较直观，可以大体看出每个功能区的主要重金属污染物；（2）运用的SPSS软件，进行主成分分析法，将多个指标转换成几个较少的互不相关的综合指标，这样使得进一步研究变得简单，避免了逐个分析这些指标而造成对研究对象的片面认识。62 模型缺点（1）在建模过程中，数据的处理比较繁琐，由于数据点多，使得运行程序时用的时间较长；（2）在绘图时没有考虑采样点的海拔，得出的空间分布图跟实际各点的分布情况不是很吻合；（3）由于所提供的金属浓度数据对应的坐标值较集中且难有内在关系，模拟的空间曲线函数与实际情况有较大的误差，而且通过预测得到的金属污染源位置有可能在实际空间中不存在，并且不能直接得出预测的金属污染源位置是属于哪个功能区。七 参考文献1 刘卫国，MATLAB程序设计与应用M，北京:高等教育出版社，2002。2 姜启源，谢金星.数学建模案例选集M，北京：高等教育出版社，2007.6。3 薛薇，统计分析与SPSS的应用（第三版）M，北京：中国人民大学出版社，2011。4 吕振通，张凌云，SPSS统计分析与应用 M，北京：机械工业出版社，2009。5 国家环境保护总局，HJ/T 166-2004土壤环境监测技术规范S，北京:中国环境科学出版社，2005。6 李小虎，汤中立，金昌市铜镍矿区周围土壤中重金属的迁移特征J，吉林大学学报：地球科学版，2009。7 张颖慧，基于环境地球化学基线的合肥地区土壤重金属元素的空间分布及污染评价，合肥：合肥工业大学，2007。8 胡克林，张凤荣，北京市大兴区土壤重金属含量的空间分布特征J，环境科学学报，2004，24(3)：463-468。附录1、问题一的程序z1=b(:,3);xi,yi,zi=griddata(x,y,z1,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(2)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(8种主要重金属元素在该城区的空间分布);text(-20,25000,570,蓝色-生活区);text(-20,25000,530,绿色-工业区);text(-20,25000,490,青色-山区);text(-20,25000,450,黄色-交通区);text(-20,25000,410,红色-公园绿地区);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,1);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(1)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(As元素浓度的空间分布);text(-20,25000,54,蓝色-生活区);text(-20,25000,50,绿色-工业区);text(-20,25000,46,青色-山区);text(-20,25000,42,黄色-交通区);text(-20,25000,38,红色-公园绿地区);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,2);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(3)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(Cd元素浓度的空间分布);text(-20,25000,2850,蓝色-生活区);text(-20,25000,2650,绿色-工业区);text(-20,25000,2450,青色-山区);text(-20,25000,2250,黄色-交通区);text(-20,25000,2050,红色-公园绿地区);z1=b(:,3);xi,yi,zi=griddata(x,y,z1,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,3);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(4)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(Cr元素浓度的空间分布);text(-20,25000,1350,蓝色-生活区);text(-20,25000,1250,绿色-工业区);text(-20,25000,1150,青色-山区);text(-20,25000,1050,黄色-交通区);text(-20,25000,950,红色-公园绿地区);z1=b(:,3);xi,yi,zi=griddata(x,y,z1,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,4);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(5)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(Cu元素浓度的空间分布);text(-20,25000,4200,蓝色-生活区);text(-20,25000,3850,绿色-工业区);text(-20,25000,3500,青色-山区);text(-20,25000,3200,黄色-交通区);text(-20,25000,2900,红色-公园绿地区);z1=b(:,3);xi,yi,zi=griddata(x,y,z1,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,5);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(6)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(Hg元素浓度的空间分布);text(-20,25000,28000,蓝色-生活区);text(-20,25000,26000,绿色-工业区);text(-20,25000,24000,青色-山区);text(-20,25000,22000,黄色-交通区);text(-20,25000,20000,红色-公园绿地区);z1=b(:,3);xi,yi,zi=griddata(x,y,z1,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,6);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(7)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(Ni元素浓度的空间分布);text(-20,25000,200,蓝色-生活区);text(-20,25000,185,绿色-工业区);text(-20,25000,170,青色-山区);text(-20,25000,155,黄色-交通区);text(-20,25000,140,红色-公园绿地区);z1=b(:,3);%高度xi,yi,zi=griddata(x,y,z1,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,7);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);figure(8)mesh(xi,yi,zi,ci)caxis(min(c4) max(c4)colorbartitle(Pb元素浓度的空间分布);text(-20,25000,650,蓝色-生活区);text(-20,25000,610,绿色-工业区);text(-20,25000,570,青色-山区);text(-20,25000,530,黄色-交通区);text(-20,25000,490,红色-公园绿地区);z1=b(:,3);%高度xi,yi,zi=griddata(x,y,z1,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);xi,yi,ci=griddata(x,y,c4,linspace(min(x),max(x),319),linspace(min(y),max(y),319),cubic);a=xlsread(D:file1.xls);b=xlsread(D:file2.xls);format short ex=b(:,1);y=b(:,2);z=a(:,8);%浓度c4=b(:,4);xi,yi,zi=griddata(x,y,z,linspace(min(x),max(x),319),linspace(min(y