高密度电阻率法25维反演演示程序使用说明c

1、高密度电阻率法2.5维反演演示程序hdrminv_ys_c使用说明书中国地质大学电法科研组 2005年7月 hdrminv_ys是用以说明高密度电阻率法2.5维反演程序使用方法和反演效果的演示程序。它使用Fortran语言设计，在QuickWin工程环境下编译为执行程序，可以在任何IBM PC或其兼容机上运行。一、输入文件 反演所需的实际观测条件参数和实测视电阻率数据，采用数据文件输入。演示程序hdrminv_ys，只能用来对给定的三个输入数据文件demo1.dat, demo2.dat和demo3.dat做反演。其中，demo1.dat是一个对称四极（施伦贝格尔）装置的正演计算理论模型数据；

2、demo2.dat是浙江地球物理研究所采用b装置测得的实测数据；而demo3.dat是物化探研究所在秦始皇陵探测中，用双向单极梯度装置的实测结果。虽然我们的高密度电阻率法2.5维反演程序，可以用来反演以下九种常用电极装置的测量数据：装置类型1¾¾双向单极梯度装置NMAMN，MN = a ，AM = n·a，B®µ。其中，a¾相邻电极间距，即测点距；n¾电极隔离系数。装置类型2¾¾正向单极梯度装置（三极装置）AMN，MN = a ，AM = n·a，B®µ。装置类型3¾

3、;¾反向单极梯度装置（三极装置）NMB，MN = a ， MB = n·a，A®µ。装置类型4¾¾对称四极(施伦贝格尔)装置AMNB，MN = a，AM = NB = n·a。装置类型5¾¾偶极-偶极装置ABMN，AB=MN= a，BM = n·a。装置类型6¾¾ 电位装置（二极装置）AM，AM= n·a， B®µ和N®µ。装置类型7¾¾ a装置（温纳装置）AMNB，AM=MN=NB= n·a。装

4、置类型8¾¾ b装置ABMN，AB=BM=MN= n·a。装置类型9¾¾ g装置AMBN，AM=MB=BN= n·a。但这里提供的演示程序只能反演上述三个数据。这些输入数据文件的共同格式是：LINENarray,Nd,AES,Nmin,Nmax H(i), i = 1,NdIs0,dIs,Xs0 Ps0其中：1） 第一行LINE为说明测区和测线信息的字符串，字符数不超过60个。2） 第二行输入观测装置的有关参数。其中： 整型数Narray为电极装置类型号。根据所用电极装置类型，取值为19之间。 整型数Nd取值为被反演的实测剖面上布置的

5、电极总数。约定剖面上第一个电极的位置序号为1，以后依此类推，直到最后一个电极的位置序号为Nd。 实型数AES =a取值为观测剖面上相邻电极的间距，也等于测点距，单位米。整型数 Nmin和Nmax分别为电极装置的最小和最大间隔系数，所谓间隔系数是指供电电极和相邻测量电极间距与a的比值n。3） 第三行输入布置电极的各测点的高程（实型数组），共Nd个，单位为米。4） 第四行输入三个整型或实型数：Is0¾第1个电极（测点）的编号，整型数；dIs¾相邻电极（测点）的编号增量，整型数和Xs0¾第1个电极（测点）的X坐标（米），实型数。 5） 从第五行开始以数组形式输入实测视电

6、阻率 Ps0。在第一个演示数据demo1.dat中，Ps0按重庆地质仪器厂所生产仪器的观测数据格式排列；在第二个演示数据demo2.dat中，Ps0按重庆奔腾仪器厂所生产仪器的观测数据格式排列；而在第三个演示数据demo3.dat中，Ps0按物化探研究所的观测数据格式排列。二、程序的运行双击执行程序hdrminv_ys.exe便可开始运行本演示程序。运行时，屏幕出现类似WINDOWS的QuickWin用户界面。为便于观察，用户需将其中的“高密度电阻率法2.5维成像系统”窗口最大化；然后，按其上屏幕提问输入相应数据：1） 屏幕显示“请键入输入文件名”¾¾用户键入demo1.d

7、at或demo2.dat或demo3.dat，然后按“Enter”键 本文后面提及“键入×××，然后按Enter键”时，将只写“键入×××”，而省略其后的“然后按Enter键”。2）“需要对实测数据做偶然误差处理吗？（Y/N）”。本程序设计了对实测视电阻率数据做“压制过失误差”和“压制偶然误差”的预处理子段。前者压制单点突变性质的过失误差，效果甚佳，且不会影响有用异常，故列为必做内容。后者采用“移动窗口光滑”算法压制具随机性质的偶然误差，有一定效果，但可能会降低对有用异常的分辨率，故列为供用户选择的处理内容。用户可根据自己的意愿，选用

8、（键入Y或y）或不选用（键入N或n）。3）屏幕显示“请键入下面内容，这些内容将写入输出文件中：” 1屏幕显示“完成计算工作的单位:” ¾¾用户键入相应字符，然后按“Enter” 键；也可不输入任何信息，直接按“Enter”键。 2屏幕显示“计算者:” ¾¾用户键入相应字符，然后按“Enter” 键；也可不输入任何信息，直接按“Enter”键。 3屏幕显示“计算机型号:” ¾¾用户键入相应字符，然后按“Enter” 键；也可不输入任何信息，直接按“Enter”键。4）屏幕显示“正在进行计算，请稍候”。5）屏幕将显示反演计算过程，依次显示

9、: 1. “实测视电阻率个数Nps=*, 反演真电阻率个数Nsigma=*, Nsigma1=*， Nsigma2=*， 背景电阻率Bkresis=*”2. “计算雅可比矩阵 正在计算，请稍候.”3. “计算雅可比矩阵的奇异值分解 正在计算，请稍候.”4. “计算雅可比矩阵的改进广义逆矩阵 正在计算，请稍候.”5. “第1 次迭代 阻尼系数=0.3 拟合方差= 0.843E-01”“第2 次迭代 阻尼系数=0.3 拟合方差= 0.499E-01”“第3 次迭代 阻尼系数=0.3 拟合方差= 0.369E-01”6）当迭代达到预定的要求时停止迭代，并在屏幕上显示当前反演结果（包括电阻率地电断面草

10、图）；然后，屏幕显示： 1“拟和方差已达到预定要求”（当拟和方差asd £ 0.01时），或 “迭代次数已达到预定要求”（当迭代次数ite = 15时），或“拟和方差在增大或改善很小”（当连续三次拟和方差的相对变化dasd £ 3%时）。2询问“对反演结果满意了吗？（Y/N）”（用户可根据屏幕显示的反演结果，判断对反演结果是否已满意了？若满意了，键入Y或y，则反演结束；若对反演结果不满意，键入N或n，则往下继续）。3若“拟和方差在增大或改善很小”，则进一步询问“需要重算雅可比矩阵吗？（Y/N）”。若键入Y或y，则返回前述第4）项（通常，只在最初几次迭代后便出现“拟和方差在增

11、大或改善很小”时，才采用“重算一次雅可比矩阵” 的方法改善反演结果）；若键入N或n，则进入下一项。4提请输入参数值：“请键入最大迭代次数”（用户按自己的意愿键入一个适当大的整数，例如，20）。5继续提请输入参数值：“请键入所要求的迭代方差” （用户按自己的意愿键入一个适当小的小数，例如，0.01）。此后，程序继续往下迭代，直到用户对反演结果满意为止。6用户对反演结果满意时，键入Y或y，屏幕显示“运算结束”。此时，用户可根据自己的需要，按WINDOWS的常规操作方法，将屏幕显示的反演结果（包括电阻率地电断面草图）保存为位图文件，或拷贝、复制到别的文件中（见下图）。三、输出文件 反演结束后，除以上

12、述位图形式给出反演结果（包括电阻率地电断面草图）外，还给出三个输出文件。输出文件名为“O1+输入文件名”的第一个输出文件是说明被反演数据、参数、反演过程和反演结果的文挡性文件。输出文件名为“O2+输入文件名”的第二个输出文件是供绘制（真电阻率）地电断面图的数据文件。输出文件名为“O3+输入文件名”的第三个输出文件，输出数据予处理（过失误差和和偶然误差校正）前后的视电阻率数据，用以观察数据予处理的效果。下面给出这三个输出文件的例子。第一个输出文件的例子（O1demo1.dat）： 高密度电阻率法反演结果文件 输入文件名: demo1.dat 测区、测线: 高密度电阻率法理论（正演）计算数据 电极

13、装置类型: 对称四极(施伦贝格尔)装置AMNB 测线上电极总数: 52 测量电极距: 5.00 最小间隔系数: 1 最大间隔系数: 17 计算单位: cug 计算者: lyz 计算机型号: p1300 计算开始时间: 20050528 162602.258 实测视电阻率个数: 561; 反演真电阻率个数: 459; 背景电阻率: 9.72第 0 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.438E+00第 1 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.637E-01第 2 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.437E-01第 3 次迭代 阻尼系数= 0.

14、30E+00 拟合方差= 0.322E-01第 4 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.248E-01第 5 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.202E-01第 6 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.173E-01第 7 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.153E-01第 8 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.139E-01第 9 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.128E-01第10 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.120E-01第11 次迭代 阻尼系数

15、= 0.30E+00 拟合方差= 0.113E-01第12 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.108E-01第13 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.104E-01第14 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.101E-01第15 次迭代 阻尼系数= 0.30E+00 拟合方差= 0.977E-02 地电断面的真电阻率（m） 9.55 14.98 13.03 12.02 11.10 10.55 10.21 10.00 9.87 9.79 10.14 11.16 11.08 11.47 10.43 10.13 9.94 9.83 9.80

16、9.58 10.48 10.35 10.43 10.21 9.98 9.84 9.75 9.87 9.72 10.86 9.98 9.70 10.23 9.81 9.72 9.67 9.86 9.79 10.14 9.42 9.28 9.75 9.91 9.60 9.58 9.78 10.41 9.84 9.83 8.93 9.28 9.79 9.47 9.48 9.72 11.17 9.59 8.93 9.12 9.06 9.58 9.31 9.36 9.72 11.35 10.02 8.07 8.28 8.68 9.20 9.45 9.23 9.77 11.19 9.30 7.94 8.1

17、4 8.41 9.07 9.23 9.09 9.78 10.99 9.36 7.78 8.06 8.49 8.95 9.08 8.94 9.83 10.92 9.11 8.15 8.11 8.30 9.20 9.04 8.77 9.95 10.37 9.17 9.08 8.29 8.44 9.72 9.36 7.90 10.05 9.81 9.12 9.48 8.54 9.47 10.75 9.41 7.88 10.12 9.56 9.67 9.25 10.58 11.08 11.94 10.08 7.63 10.22 9.37 9.30 11.21 13.47 13.39 12.38 9.1

18、9 7.38 10.36 9.02 12.27 13.62 14.76 15.90 10.51 7.51 6.95 10.52 8.97 10.90 11.74 11.80 9.57 8.40 6.28 6.38 10.49 10.30 7.82 5.14 6.50 5.57 5.09 5.82 6.17 10.01 11.53 2.17 1.58 2.84 3.03 3.66 4.89 6.33 9.34 8.19 0.91 0.83 1.47 1.93 2.88 4.61 6.13 9.23 9.17 0.87 0.71 1.16 1.74 3.00 4.67 6.51 9.90 9.33

19、 0.87 1.02 1.49 2.20 3.58 5.12 6.99 10.59 9.58 1.75 2.11 2.48 3.03 4.12 6.20 7.78 10.85 11.60 3.36 4.02 3.67 4.23 5.17 6.97 8.45 10.38 14.22 6.58 6.11 5.40 5.72 6.52 7.96 8.93 9.73 12.98 9.26 8.24 7.06 7.36 8.17 8.76 9.37 9.57 10.19 11.32 10.10 9.36 9.40 9.45 9.53 9.24 9.70 7.70 14.67 11.68 11.62 11

20、.38 10.40 9.71 9.39 9.30 9.66 20.47 15.55 15.03 13.39 11.38 9.90 9.68 8.50 15.90 23.16 24.36 19.51 15.04 12.06 10.27 9.72 8.00 22.15 42.03 38.82 23.97 17.60 12.94 10.42 9.95 8.32 22.90 56.83 47.55 30.90 20.04 13.82 10.87 10.22 9.74 14.79 49.90 50.08 35.85 21.73 15.26 11.51 10.55 11.94 9.21 31.07 39.

21、15 33.27 23.17 15.68 11.88 10.77 13.21 6.00 17.03 25.19 25.87 20.18 15.82 12.80 11.15 12.16 6.85 12.17 15.76 15.96 15.96 15.04 13.37 11.42 10.65 8.43 11.02 12.35 11.15 11.41 13.67 12.87 11.64 9.82 10.19 10.27 10.95 10.39 9.99 11.16 12.07 11.52 9.55 11.01 9.58 11.62 10.89 9.64 9.96 11.21 11.03 9.62 1

22、0.83 8.63 12.44 12.41 8.99 9.12 10.55 11.01 9.79 10.43 7.84 12.05 12.99 9.32 8.89 10.02 10.86 9.97 10.12 8.33 9.48 12.60 10.15 8.92 9.74 10.89 10.09 9.96 8.57 8.29 12.04 10.80 9.19 9.96 10.91 10.12 9.99 8.76 8.11 11.39 10.94 9.31 9.98 10.95 10.08 10.19 9.34 8.25 10.75 11.04 9.82 10.04 10.56 10.03 10

23、.29 9.85 8.73 10.66 10.91 10.12 10.54 10.21 10.00 10.33 9.91 9.32 10.57 10.94 10.39 10.59 10.24 9.99 9.90 10.64 9.60 10.55 11.00 10.34 10.57 10.26 10.00 9.80 10.87 10.67 10.47 11.06 10.74 10.39 10.24 9.98 9.83 11.19 10.97 10.82 10.68 10.50 10.34 10.21 9.75 10.65 11.58 11.02 10.81 10.60 10.43 10.28 1

24、0.16计算开始时间： 20050528 162602.258计算结束时间： 20050528 162733.619上述第一个输出文件例子中，前13行表示被反演资料的条件参数，第14行输出反演的有关参数，供参考。第15行开始显示迭代过程¾¾各次迭代的阻尼系数和拟和方差。此后，输出反演所得地下岩石的真电阻率¾¾第一行数据表示地电断面从左算起第一列的电阻率，其余依此类推。这种输出形式便于有经验的用户一目了然地看出反演所得地下电性分布。在输出文件的最后两行，分别输出了开始反演和结束反演的时间，由此易于算出反演计算花费的时间。第二个输出文件的例子（O2demo1

25、.dat）： 测点号 x（m） z（m） 真电阻率(m) 1.50 2.50 -2.50 9.55 1.50 2.50 -7.50 14.98 1.50 2.50 -12.50 13.03 1.50 2.50 -17.50 12.02 1.50 2.50 -22.50 11.10 1.50 2.50 -27.50 10.55 1.50 2.50 -32.50 10.21 1.50 2.50 -37.50 10.00 1.50 2.50 -42.50 9.87 2.50 7.50 -3.00 9.79 2.50 7.50 -8.00 10.14 2.50 7.50 -13.00 11.16 2.

26、50 7.50 -18.00 11.08 2.50 7.50 -23.00 11.47 2.50 7.50 -28.00 10.43 2.50 7.50 -33.00 10.13 2.50 7.50 -38.00 9.94 2.50 7.50 -43.00 9.83 3.50 12.50 -4.00 9.80 3.50 12.50 -9.00 9.58 3.50 12.50 -14.00 10.48 3.50 12.50 -19.00 10.35 49.50 242.50 -12.50 10.47 49.50 242.50 -17.50 11.06 49.50 242.50 -22.50 10

27、.74 49.50 242.50 -27.50 10.39 49.50 242.50 -32.50 10.24 50.50 247.50 7.50 9.98 50.50 247.50 2.50 9.83 50.50 247.50 -2.50 11.19 50.50 247.50 -7.50 10.97 50.50 247.50 -12.50 10.82 50.50 247.50 -17.50 10.68 50.50 247.50 -22.50 10.50 50.50 247.50 -27.50 10.34 50.50 247.50 -32.50 10.21 51.50 252.50 7.50

28、9.75 51.50 252.50 2.50 10.65 51.50 252.50 -2.50 11.58 51.50 252.50 -7.50 11.02 51.50 252.50 -12.50 10.81 51.50 252.50 -17.50 10.60 51.50 252.50 -22.50 10.43 51.50 252.50 -27.50 10.28 51.50 252.50 -32.50 10.16 1.50 2.50 0.00 99999.00 2.50 7.50 -0.50 99999.00 3.50 12.50 -1.50 99999.00 4.50 17.50 -2.50

29、 99999.00 5.50 22.50 -2.50 99999.00 6.50 27.50 -1.50 99999.00 7.50 32.50 -0.50 99999.00 8.50 37.50 0.50 99999.00 9.50 42.50 1.50 99999.00 10.50 47.50 3.00 99999.00 11.50 52.50 5.50 99999.00 12.50 57.50 8.00 99999.00 13.50 62.50 10.00 99999.00 14.50 67.50 11.50 99999.00 15.50 72.50 12.50 99999.00 16.

30、50 77.50 12.50 99999.00 17.50 82.50 11.00 99999.00 18.50 87.50 9.00 99999.00 19.50 92.50 7.00 99999.00 20.50 97.50 4.50 99999.00 21.50 102.50 1.50 99999.00 22.50 107.50 -1.50 99999.00 23.50 112.50 -4.50 99999.00 24.50 117.50 -7.50 99999.00 25.50 122.50 -10.00 99999.00 26.50 127.50 -12.00 99999.00 27

31、.50 132.50 -14.00 99999.00 28.50 137.50 -15.50 99999.00 29.50 142.50 -16.50 99999.00 30.50 147.50 -16.50 99999.00 31.50 152.50 -15.50 99999.00 32.50 157.50 -14.00 99999.00 33.50 162.50 -12.00 99999.00 34.50 167.50 -10.00 99999.00 35.50 172.50 -7.50 99999.00 36.50 177.50 -4.50 99999.00 37.50 182.50 -

32、2.00 99999.00 38.50 187.50 -0.50 99999.00 39.50 192.50 0.50 99999.00 40.50 197.50 1.50 99999.00 41.50 202.50 2.50 99999.00 42.50 207.50 3.50 99999.00 43.50 212.50 4.50 99999.00 44.50 217.50 5.50 99999.00 45.50 222.50 6.50 99999.00 46.50 227.50 7.50 99999.00 47.50 232.50 8.50 99999.00 48.50 237.50 9.

33、50 99999.00 49.50 242.50 10.00 99999.00 50.50 247.50 10.00 99999.0051.50 252.50 10.00 99999.00第二个输出文件逐点输出地电断面中各点正上方的测点编号Ix,横向（X）坐标，纵向（Z）坐标和电阻率值。这样的输出数据，可直接供Surfer等绘图软件绘制地电断面（电阻率）等值线图。其中，最后Nd-1=51行数据，是专为绘制不平地面线用的。第三个输出文件的例子（O3demo1.dat）： 过失误差和偶然误差校正前后的视电阻率数据对比 输入序号I A极序号 极距序号 校正前Ps(I) 校正后Psc(I) 1 1 1 9.75 9.67 2 2 1 9.63 9.87 3 3 1 10.04 9.89 4 4 1 10.00 9.79 5 5 1 9.33 9.73 6 6 1 9.70 9.53 7 7 1 9.53 9.44 8 8 1 9.