张玉君
（地矿部航空物探遥感中心）
摘要：本文概述伽玛能谱测井数据处理、解释、图形图像软件系统功能、设计思想和研究特点、环境和安装、及系统主要测试结果。提出“雏形程序”的思想，雏形程序是细部设计阶段的得力工具。
序言
本系统用C语言编写，可在IBM PC/286及以上机型或兼容机上运行，具模块化结构。研制过程中严格以软件工程方法为指导，历经总体设计、细部设计、编程和调试、系统集成和系统测试等阶段。除原设计中预期完成的正反演定量解释、趋势分析、彩色图形显示和交会图3组20个功能外，还扩展增补了预处理和交互解释两组8个功能，特别是剖面对比交互解释、岩性柱交互解释、综合显示和岩性柱交互解释功能的产生使整体水平升高，形成更为完美的交互解释系统。
本系统旨在解决自然伽玛能谱测井定性定量解释的计算化问题。我部对普查急缺矿种钾盐的方法技术研究给予极大支持，在20世纪60年代物探所就提出了用自然伽玛能谱测井找钾的方法，经过方法试验和五个型号仪器的研制，在六个省得到推广应用，取得一定生产效果，本系统的研制成功无疑使我国钾盐能谱测井技术得以完善化。
一、系统功能
伽玛能谱测井数据处理、解释、图形图像软件系统（HRSIS）包括28个功能，其主要用途按组依次分述如下：
A:PRP（A组：预处理）
A.1 CALCALT（计算），对两条测井曲线进行算术运算。
A.2 CALCULT1（计算1），对两条测井曲线进行算术运算，输出格式与输入格式相同。
A.3 LGMINMAX（测井曲线最大最小值），查找测井曲线的最大最小值。
A.4 LITHOLOG（岩性柱），形成岩性柱数据文件，用于绘制测井剖面图。
A.5 SEPARATR（曲线分离），将多条曲线的数据文件分离为多个单条曲线的文件，其格式适于进一步处理或显示。
B:KTU（B组：K、Th、U含量计算）
B.1 W2F（双窗法求含量），利用双窗法计算测井曲线各点钾元素，干扰元素（U+Th）及其他元素的含量。
B.2 W2FF（双窗法求K含景），利用双窗法计算测井曲线各点钾元素的含量。
B.3 W2_GJ（解方程组求双窗法系数），利用GaussJordan函数求解实型线性方程组，得出模型标定的双窗法系数。
B.4 W2_MR（多重回归拟合求双窗法系数），利用MultipleReg函数对观测数据集建立一个最小二乘线性多重回归方程，解出模型标定的双窗法系数。
B.5 W3F（三窗法求含量），利用二窗法计算测井曲线各点钾元素、钍元素、铀元素及其他元素的含量。
B.6 W3_GJ（解方程组求三窗法系数），利用Gauss-Jordan法求解线性方程组，得出模型标定的三窗法系数。
B.7 W3_MR（多重回归拟合求三窗法系数），利用MultipleReg函数，对观测数据集建立最小二乘线性多重回归方程，解出模型标定的双窗法系数。
C:TRD（C组：趋势分析）
C.1 TRENDML（多重趋势分析），利用DataSmoothWeights函数及输入的加权系数对测井曲线进行多重趋势分析。
C.2 TRENDSG（GS趋势分析），利用DataSmoothSg函数进行趋势分析，该函数采用的算法是简化最小二乘平滑和微分法的Savitzky-Golay技术。平滑方式有5点、7点、9点、11点和13点五种；对于9点平滑方式并计算一阶导数时，程序做特殊处理，采用待定系数。
C.3 TRNDSG1（GS趋势分析1），利用Savitzky-Golay技术对测井曲线进行趋势分析。
C.4 TRENDW（加权趋势分析），利用DataSinoothWeights函数及输入的加权值对测井曲线进行趋势分析。
D:GRP（D组：图形显示）
D.11LEGENDS（图例1），在页面下部创建14个岩性彩色花纹图例。
D.22LEGENDS（图例2），创建14个岩性彩色花纹图例。
D.3 KTHUCTBW（K、Th、U含量黑白图），绘制测井钾、钍、铀元素含量黑白图。
D.4 KTHUCTCL（K、Th、U含量彩色图），绘制测井钾、钍、铀元素含量彩色图。
D.5 KTUCTBW1（K、Th、U含量黑白图），同D.3，但Th和U曲线从右边起始。
D.6 LOGMTPBW（多条测井曲线黑白图），绘制多达10条测井曲线黑白图。
D.7 LOGMTPCL（多条测井曲线彩色图），绘制多达10条测井曲线彩色图。
D.8 SCATTGRM（黑白散布图），利用库函数＿pg_chartscatter对两条测井曲线做黑白散布图。
D.9 SCATTG15（彩色散布图）对两条测井曲线做彩色散布图。
E:ITR（E组：交互解释）
E.1 CORRELAT（对比交互解释）。交互绘制多钻孔测井柱对比图。
E.2 LOGMTP（岩性柱交互解释），显示多条测井曲线，并对岩性柱进行交互解释。
E.3 VICTORY（综合显示及岩性柱交互解释），显示多条测井曲线、3条含量曲线、2条趋势曲线，并对岩性柱进行交互解释。
二、系统设计思想及研制特点
1.以软件工程理论为指导
本系统研制自始至终注意和强调以软件工程理论为指导，提高软件研制水平。经历了立项阶段，总体设计阶段、细部设计阶段、编程和调试阶段、以及总体测试和总结阶段。现将通常最易忽视的前三个阶段简述一下。
在立项阶段进行了可行性研究和系统功能设计，由于我们对方法有较深入的了解，在此阶段将系统所涉及的教学方法和公式提得足够明确，这阶段的主要成果是课题设计书。
在总体设计阶段进行了结构设计、用户接口设计、以及针对所需计算方法选择适用的软件资源。在结构设计中从课题设计书中所定系统总体功能出发，用自顶向下、逐层分解和模块化的设计方法，设计了系统的层次和构成模块，系统结构在软件的开发过程中又经自底向上补充和修改，形成最终的层次体系。此阶段还从丰富的软件资源中选择了可借用函数，如：GaussJordan解方程组函数、MultiplReg多重回归拟合函数、DataSmoothWeights和DataSmoothSg曲线光滑函数、＿рg_chartsca-ter两维概率统计图形库函数等。
在细部设计阶段主要进行模块设计，它的任务是在原理设计和结构设计的基础上，针对每个模块设计实现其指定功能、算法和外部接口要求的内部数据结构和程序型逻辑结构，编写和调试“雏形程序”。“雏形程序”的提出是在本系统研制过程中我们自己创造的一种做法，针对每一类功能近似的模块，设计一种雏形程序，它在逻辑上是完善的，可运行的，但在细节上是有待补充的。我们体会到雏形程序对于编程好处极大，一个雏形程序一般仅数十条程序语句，而一个完善的最终主函数（程序）可以长达700多条语句并调用十多个函数（子程序）；后者以前者为基础，从简到繁循序渐进逐步增补完善，简化了复杂程序的调试，节省了时间，提高了源程序文本的质量。大多数雏形程序都保存下来，留做系统维护时参考。
2.模块化软件结构
HRSIS系统具有层次性和模块化的软件结构。系统深度为4层，用户可见的有3层（见图1），第一层为总控模块。第二层功能组合，第三层为28个主要功能模块.第四层为公用模块。整个系统由98个C语言编写的函数组成，其中主函数有29个。全都源程序共有语句一万多条共占1.165M字节，正好容纳在一片l.2M的软盘上。总控模块（HRSIS）给出系统入口界面，并通过spawn函数创建并执行子进程，每一个功能模块为一个独立的子进程，表现为一个.EXE文件，可被spawn调用。而HLP信息则通过.BA T批处理文件调用，因此为每一个HLP文件提供一个.BAT文件，每一个批处理也当做一个子进程被spawn函数调用。
图1 HRSIS系统用户界面
3.友好的用户接口
HRSIS系统具有友好的用户接口，它有以下特点：
*使用简便，全部功能由一个总控模块所控制，进入入口程序后，选用功能时，用户只需键入序号即可，无需重复功能名。
*请用户键入参数的提示术语是标准的一致的，并在括号内均给出示范。
*交互功能均有操作提示。
*拥有完善的HLP子系统，各功能的用法及操作实例都可从 HLP子系统获得。此外为了系统维护（纠错性维护、扩充性维护更新）的方便，在HLP子系统中还包括有多数子程序（函数）的HLP文件。为了英语困难的用户还提供有中文HLP子系统，它需在中文DOS下运行。
4.采用分段、分页处理，解决宏数据量与内存不足的矛盾
本系统无论计算性功能或是图形显示性功能均面临大规模数据的处理，这些数据都存在数据文件内，每一条测井曲线数据文件可以包括数千实型数据点，对于图形显示功能往往要同时显示10多条曲线，由于内存的限制，不可能一次把多条曲线都调入内存，为了解决宏数据量与内存不足的矛盾，我们采用了动态分配数值和分段，分页数据处理技术。
5.建立完善的规范化的文档
为了提高软件产品的可维护性和可继承性，本系统在研制过程中十分重视文档的建立和完善。各阶段都有文字材料。特别是进入编程阶段之前，尚在细部设计阶段就制定了源程序的编写规则，特别强调注释的运用。使不同人员不同时间所编写的源程序都保持统一的程序风格，且文本美观、清晰、易于理解。中、英文的HLP文本都是在各函数通过后立即按统一格式编写，这就大大简化了后期工作。
本系统提交的文档计有：
*“伽玛能谱测井数据处理、解释、图形图像软件系统”课题没计书。
*HRSIS系统使用手册。
*HRSIS系统源程序集。
*HRSIS系统源程序安装软盘。
*HRSIS系统研制报告。
三、系统环境及安装
1.系统硬、软件环境
主机：该系统在IBM PC 286上开发成功，并在IBM PC 386上安装。主存1M以上，显示器为640×480高分辨率彩显，带VGA视频卡和1.2M软盘驱动器。
打印机：EPSON LQ-1600K。
操作系统：MS-DOS 3.31以上版本，中文DOS。
语言工具：Quick C。
2.系统安装
HRSIS系统全部源程序共占1.165MB，容纳在一张1.2MB软盘上，它下设11个子目录：HR-SIS、HLP、BAT、PRP、KTU、TRD、GRP、ITR、FNCT、SMPSDL和SMPZYJ。此外还有一个启动程序：SUHRSIS.BAT。
系统安装时，首先置于系统未来所处的目录下，运行下述两条命令即可建立系统目录树并把系统源程序装入所建立的目录树中。
CORY A:SUHRSIS.BAT←
SUHRSIS←
进行编译连接后，使可通过HRSIS.EXE总控函数运行HRSIS系统了。
四、系统测试
我们对软件测试基本上分两个步骤进行的，即模块测试（或称单元测试）和系统测试（或称集成测试）。
我们在单元测试努力中采用穷尽法，即对模块的每一逻辑结构、各种可能的输入、以及各种可选的操作选项进行详尽的测试。之后又在IBM 286和IBM 386机上进行了四次系统集成，两套平行集成测试，系统可靠性得到充分证明。现将系统测试中有代表意义的结果列举如下。
1.定量解释功能
定量解释（KTU组）共有七个功能模块，其中四个为求定量解释系数，即正演；三个为根据所求之定量换算系数计算放射性元素含量，即反演。根据放射性测井和钾盐能谱测井的生产实际所用三窗法和双窗法两种方法，本软件系统具有三窗法和双窗法的正反演功能。
三窗法的关系式在设计书中有详细说明，以矩阵形式可表示为；
N=AQ
式中：计数率向量N为经本底修正了的K、U、Th窗的计数率；
含量向量Q为K、U、Th的岩层含量；
A为灵敏度矩阵。
正演的任务利用标定数据N和标定模型已知含量Q，求解灵敏度矩阵A。
反演的任务是根据已知的A矩阵和计数率N，求含量向量Q。
双窗法是用于钾盐能谱测井的一种方法，由于U、Th元素含量甚低，把它们合并视为干扰元素。双窗法关系式虽不同（见设计书），但以矩阵形式表达同上式。
在系统设计中为解此式采用了两种数学方法：最小二乘线性多重回归拟合法（MultipleReg）和解线性方程组的GaussJordan法。针对双窗法和三窗法，形成四个正演模块：W2_MR、W3一MR、W2_GJ和W3_GJ。从数学原理上讲后一种方法更为严密，但对实测数据却是第一种方法更实用、更准确。我们建议一般均应使用：W2_MR和W3_MR模块，只有当模型数量少（对双窗法最少为二个，对三窗法最少为三个）的情况下才不得已使用W2_GJ和W3_GJ模块，为此把这两个功能保留下来了。
W2F、W3F和W2FF为反演模块，我们选用了1988年在石家庄大郭村静态板上对多道航空能谱仪的实际标定数据作为测试数据。
表1 大郭村机场航空能谱静态模型标准含量
列入验证数据表中的数据为两箱探测晶体经过死时间校正扣除仪器本底的净峰面积，对双窗
法，我们将U、Th窗计数相加作为N干，模型板Qu,QTh相加做为Q干。
表2 三窗法验证数据
表3 二窗法验证数据
此外我们还准备了测试含量计算用的测井数据文件。
用表2的数据运行W3_MR功能，所得A阵为：
113.63-0.87-0.59
-5.365.07-0.50
1.33-0.505.12
再用灵敏度系数文件模拟测井数据文件运行 W3F，得到 AP-M标定板 K、U、Th含量为：Qk=3.93%；Qu=20.59ppm；QTh=40.50ppm。
用表3数据运行W2_MR功能，所得灵敏度系数矩阵如下：
113.72-0.63
-4.024.61
再用此灵敏度系数文件及模拟测井数据文件运行W2F，得到AP-M标定板K和干扰元素含量为：Qk=3.89%；Q干＝61.96ppm。与表1已知数据对比，结果准确，完全可以证明这些功能运算是正确的。
2.测井曲线图形显示及交互解释功能
测井曲线图形显示在GRP组的D.3-D.7五个功能中体现，并提高成交互解释功能LOGMTP和VICTORY。验证VICTORY功能对这些功能有代表性；而且VICTORY中还显示了趋势分析曲线及图例，也就同时验证了TRD组功能及D.1和D.2功能。
我们选云南省兰坪县拉井盐矿ZK4号钻孔能谱测井结里为验证数据，原始能谱测井曲线有4条：右积分曲线、左积分曲线、峰顶微分曲线和自然伽玛曲线，经数字化，趋势分析和岩性柱数据分析，运行VICTORY功能。图2（彩版附图）展示出此功能漂亮的结果，共显示九条曲线（三条含量曲线、四条实测曲线、两条趋势分析曲线）两个测井柱，其中一个可以进行交互修改。在274米和325米处显示出明显的钾盐岩层。最后页（图3（彩版附图1））下部为岩性图例。
3.对比交互解释功能
CORRELAT功能允许我们对多个钻孔部面做对比解释，最多可以显示10个钻孔部面，对比层数不受限制。屏幕下部给出了操作提示。钻孔剖面可以用含量曲线、原始曲线或是趋势曲线。
图4（彩版附图1）为验证结果。
4.两维概率统计功能
图5（彩版附图1）为D.9 SCATTG15彩色散布图功能的验证结果，不同概率的点以不同色彩和不同符号表示，本例中最高重复概率为69，这意味着有69个点的数据在两条测井曲线上分别具有接近的值。
这种图也称为交会图，是测井分析的重要方法特别是对于岩相的研究十分有用，在图5的上部的点群代表着含Th最高而含K量中、低的岩性，很可能是泥质岩。
致谢 本系统研制过程中，史殿林、谢欣、李昌国等同志参加部分工作，特此致谢。
参考文献
[1]章晔，华荣洲，石柏慎.放射性方法勘查.北京：原子能出版社，1990
[2]李朝栋：校准核测量仪器用的标准系列模型的活度参数.航测与遥感，1988（4）
[3]张海藩：软件工程导论.北京：清华大学出版社，1992
[4]崔俊芝，黄玉霞，韩其瑜.软件工程方法.北京：科学出版社，1992
[5]Microsoft C科学与工程工具库，北京希望电脑公司，1991
[6]最新Microsoft Quick C使用大全，卷1、2、3、4，中国科学院希望高级电脑技术公司，1990
[7]John H.Doveton,Log analysis of subsurface geology—concept and computer methods,USA,John Wiley&Sons.Inc.,1986.
[8]Spectralog,USA,Dresser industries,Inc.,1981.
[9]Автоматизированная интегрированная ивформационно－обрабатывающая системадля геоφизическихисследований неφтяцых и газовых скважин，Россия，1992.
GRAPHIC AND IMAGE SOFTWARE SYSTEM FOR GAMMA-SPECTRALOG DATA PROCESSING AND INTERPRETATION
Zhang Yu jun
（Center of Aero-Geophysics and Remote Sensing,MG&MlR）
This paper briefly describes the functions,the design ideas,the development specjalty,the hardware and software environment,the installation,and themost important test results.The new idea was given about the“EMBRYONlC PROGRAM”，whjch can be the important tool forthe detailed design stage.
彩图见彩版附图1