坑透数据处理

无线电波透视法（坑透）从上世纪20年度苏联A·A·彼得罗夫斯基院士提出设想至今，在坑道找矿、煤层工作面构造探测等方面取得了很好的应用效果。对于无线电波透视法的资料处理，早期采用手工计算、绘图，比较繁琐、工作量大；随着电子计算机技术的发展，采用计算机自动处理的技术逐渐被推广应用，使资料处理效率得到了很大的提高，减轻了地球物理工作者的劳动强度。从无线电波透视法理论出发，结合无线电波资料处理、解释以及生产报告撰写要求，设计了条件试验分析、异常数据的编辑、多方法的求参分析功能、单点综合曲线绘制、单巷道多测点综合曲线绘制、工作面内吸收系数的反演计算等功能，满足无线电波透视法资料处理高效、准确的要求。坑透资料的基本处理流程如下：一、参考处理流程对坑透资料的基本处理步骤如下：1.进行“文件-新建工程”操作，输入工程名称、各拐点坐标等工程参数。2.进行“数据-测值输入”操作，输入各发射点对应的实测场强值及参数。3.若存在条件试验，可以执行菜单中的“条件试验分析”，采用回归算法分析初始场强H0及吸收系数β。4.进行“参数分析”操作，可以选择“实测场强法”、“图解法”、“回归分析法”，对各个发射点对应的接收场强值进行初始场强H0及吸收系数β分析，掌握测区的电磁参数。5.进行“层析成像计算”操作，选择“绝对衰减法”或“相对衰减法”进行CT成像反演计算。6.进行“结果输出”操作，将反演结果输出为文本格式，可以利用Surfer等专业绘图软件进行绘图操作。7.进行“图形查看”操作，可以根据资料解释、报告撰写的要求，输出各种曲线。对于已经存在的工程，可以执行“打开工程”操作，然后再进行3~7步。二、新建工程点击菜单中的“文件”-“新建工程”，弹出对话框：图2-1新建工程对话框http://好pan.好baidu.好com/s/好1o6I5uyu——去掉“好”，见演示程序。（1）输入工程名称、组成工作面的拐点总数、发射点数目。（2）输入4条巷道的测点参数。巷道：对于一个工作面，规定由4条巷道组成，即上、下巷、切眼、联络巷。巷道编号按顺时针进行，上巷为1号巷。。。。。拐点：即组成工作面的各条线段的端点。拐点按顺时针编号，左上角第一个点为0号点。测点编号：每个巷道的测点编号独自进行。上、下巷从左往右编号，点距10米，起点为0号，向左1、2、3…以此类推。切眼和联络巷从上往下编号，起点为0号。图2-2巷道、拐点编号（3）输入拐点坐标。对于拐点坐标约定：(a)垂向上最上方一个点Y坐标为0；（b）横向上最左面一个点的X坐标为0。由以上两点确定坐标原点，X坐标向左为增大方向，Y方向向下为增大方向。坐标值都大于零。如图2-3所示。图2-3拐点坐标规定三、打开已有工程点击“文件”-“打开工程”菜单，可以打开已保存的工程，进行重新计算或查看、出图操作。四、数据“数据”主菜单下包括“测值输入”、“测值修改”、“工程参数修改”三个子菜单。（1）测值输入点击“测值输入”菜单后，弹出4-1对话框，输入逐个发射点对应的观测值。观测值由逗号隔开。输入值后可以“检查”输入的数据个数是否与输入接收点匹配。图4-1定点观测数据输入（2）修改数据点击“修改数据”子菜单，可以选择对已经输入的数据进行修改。图4-2为选择待修改的记录编号对话框；图4-3为修改数据对话框。图4-2选择修改记录号图4-3修改数据对话框（3）修改工程参数若工程参数输入有误，可以对其进行修改，点击“工程参数修改”子菜单，弹出图4-4所示对话框。需注意：修改工程参数后要重新输入测值。图4-4修改工程参数对话框五、编辑在数据输入后，可以对单独发射点的接收数据中的畸变数值进行编辑。点击“编辑”菜单，客户区生产单发射点对应的接收曲线，数据点用红点表示，可以通过鼠标进行数据值的修改。鼠标点击数据点，数据点被选中后变成绿点，移动鼠标到下点，点击屏幕，数据点将移到新位置，同时数据值也在内存中改变。如图5-1、5-2所示图5-1被选中的数据点图5-2修改后的数据点六、条件试验该项功能是对条件试验所观测的数据利用线性回归法来求取初始辐射场强H。和煤层的吸收系数β。点击“条件试验”菜单，弹出条件试验数据录入对话框，如图6-1所示。图6-1条件试验数据录入利用线性回归算法获得求取初始辐射场强H。和煤层的吸收系数β。H。=52.64dB，β=0.0254dB/m。并利用计算所得的初始场强和煤层吸收系数计算各个接收点的理论场强值H0。图6-2为分析结果，图中蓝线为实测值，红虚线为计算所得理论场强曲线，棕色线为衰减系数曲线。图6-2条件试验回归分析结果七、参数分析因为条件试验的数据有限，通常在无线电波透视工作时往往利用所有测量值进行分析，通过大量的实测值分析，更能找到适合工作区的电磁参数。对于电磁参数的分析一般采用三种方法：实测计算法、图解法、回归分析法。（1）实测计算法选择均匀介质地段设置发射点，距离发射点为若r1和r2两点接收场强为H1和H2，由磁场表达式可知：sin1011reHHrsin2022reHHr移项相除，并取对数)(ln122211rrrHrH则吸收系数β为：122211lnlnrrrHrH图7-1实测法计算初始辐射场强H。和煤层的吸收系数β（2）图解法由于煤层并非完全均匀，采用选两点实测场强计算出H。和β与实际仍存在大的差异，为了求取准确的H。和β，可以采用多点测值用图解法求取H。和β值。图解法的原理：将多点实测的H1，H2，….Hn值乘上相应的发-收距离r1，r2，…rn作为纵坐标，r1，r2，…rn作为横坐标，作单发射点实测场强曲线图。由公式122211lnlnrrrHrH可知，理论上均匀介质形成的lnHr-r曲线应该是条直线。我们通过多点实测散点中拾取一条与实测值比较吻合的直线，直线的斜率就是吸收系数β，而r=0时的纵轴的截距就是H0。如图7-2所示。点击“参数分析”-“图解法”，在客户窗口显示0号记录的实测值lnHr-r曲线，点击左上角工具栏的按钮，可以用鼠标在图上绘制直线，并计算H。和β值。图7-2图解法获取初始辐射场强H。和煤层的吸收系数β（3）回归分析法采用图解法计算初始辐射场强H。和煤层的吸收系数β具有很大的随意性，为了取得好效果，可以采用直线回归分析方法，采用线性回归迭代算法拟合出最佳直线，求取合理的电性参数。在干扰小、煤层内部无构造情况下，回归分析法可以取得很好的效果，但是对于构造多、干扰多的地段，回归分析法效果不好。回归分析法原理：线性回归法是用一系列测点，进行最小二乘法意义下的线性拟合。设有n个实测的实iH，i=1,2,...,n,令xriiiiirHylg20实对于均匀介质，n个实测值的lnHr-r曲线理论上是条曲线，因此对于这些散点所在的曲线可以采用直线回归方程表示：bxay对于散点的直线回归方程的求解采用最小二乘法。最小二乘法的出发点是使实际测量数据yi与拟合直线上对应的估计值的残差的平方和为最小。即为使值最小，只要使a和b的偏导数为零，即可解得a和b的值。ayxxxyxxiiniiiniin21121bxynxyxxiiiniin121式中：xnxynyiiniin1111由于实测数据受构造、干扰等很多因素控制，有的散点偏离度比较大，其参与计算求出的直线回归方程与实际情况误差较大，为了获得最佳的拟合直线，需要采用迭代法进行拟合，每次迭代剔除离散大的数据，最终参与计算的散点都是紧密在回归直线附近的，其计算所得的H。和煤层的吸收系数β比较准确。图7-3线性回归法计算初始辐射场强H。和煤层的吸收系数β八、层析成像计算上世纪70年代中期，美国已故的测井学家R.J.Lyele等人率先利用直射线理论，把医学CT引入地学领域，把透射层析应用于跨孔电磁波探测的资料处理，推动了全球范围内地下物探层析技术的应用和研究。层析成像技术在无线电波透视方法中的应用，提高了该方法解释的精度，结合传统的综合曲线分析，有效解决煤层中的构造问题。图8-1反演参数（1）绝对衰减法在无线电波透视中，由于综采工作面煤层的倾角都小于25°，发射天线与观测点方向近乎正交，即在θ＝90°时，电磁场的计算公式为：reHHr0此式中H、H0的单位为微伏，ｒ的单位为米，β的单位为奈培／米。对上式两边取对数，经变换后可得：βr＝H0－H－20lgr式中H、H0的单位为分贝，r的单位为米，β的单位为分贝／米。如图所示是一个工作面被网格化后的示意图。把网格化后每个均匀的小块称为一个像素，在此区域内有一条射线yi穿过了衰减系数分别为xxxn12,,...,的诸像素，并在这些像素上的截距分别为dddiiin12,,...,。这样在第i条射线路径上则有:rdxijjjn1把式（1－5）代入式（１－４）中可得到第ｉ条射线的方程：dxyijjijn1这里yHHriii010lg式中Hi──第i次观测的实测场强值。ri──第i条射线长度，即rdiijjn1。若进行多重观测，即用多个波源发射，对每个发射点进行多点接收测量，即可得到如下矩阵方程：dddddddddxxxyyynnmmmnmm1112121222121212也就是DX=Y式中：D──mn阶系数矩阵，为射线穿过诸像素的截距。其中D的元素dij是第i次观测中射线被第j个像素的截距。i＝1，2…m，m是射线总条数（观测次数）。j＝1，2…n，n是工作面网格化后像素的个数。X──是n1阶未知数矩阵，代表诸像素的绝对衰减值。Y──m1阶常数矩阵，是各观测方式下与实测场强有关的常数。由于方程的建立是基于各像素绝对衰减作用的叠加，因此求得的各像素的绝对衰减系数也就体现了介质的物性。对于方程组D*X=Y的求解，基本上可分为两大类，一类是基于付里叶变换或Radon逆变换的变换方法，一类是代数方法。后者又可分为矩阵反演法和迭代重建法。属于矩阵反演的算法有：SVD（奇异值分解），CCG（约束共扼梯度法），LSQR（最小二矩阵分解）；属于迭代重建法的有ART（代数重建技术），SIRT（联合迭代技术）和BPT（反向投影技术）。图8-2绝对衰减反演成果（2）相对衰减法相对衰减层析成像的控制方程：DXS（1－10）SHrrHiiii010lg（1－11）xxjj（1－12）这里β是介质背景衰减系数，D是mn阶距离系数矩阵，ΔX是n1阶未知数矩阵，S是m1阶常数矩阵。ΔX的物理意义是测区内某处的绝对衰减系数与该处背景衰减系数之差。因此ΔX可称为“相对衰减”，类似地可把式称为“相对衰减方程”，解此方程得到的是地下介质相对衰减的二维分布。相对衰减层析成像的控制方程与波源的初始辐射强度无关，仅与背景场有关。由于在复杂地电条件下也能确定背景场，所以大大扩展了电磁波层析成像的应用范围。相对衰减层析成像控制方程的建立不受波源初始辐射强度变化的影响，并且能够抑制部分数据噪声，同时重建过程自动消除了不均匀的背景衰减，所以重建图像的可靠性和分辩率较高，背景清晰，异常突出。图8-2相对衰减反演成果九、层析成像结果输出点击“成果输出”菜单，可以输出绝对衰减、相对衰减反演到文件中，以方便用Surfer等软件绘制反演成果。数据输出分三列，分别为：X坐标Y坐标吸收系数值由于有的工作面试不规则形状的，需要对其进行白化，点击“白化”，可以输出白化文件。十、图像查看在资料解释中和成果报告撰写中都需要用到综合曲线图等图件，软件中设计了单发射点实测曲线、综合曲线，上下巷实测曲线、综合曲线等图件绘制功能。（1）射线图点击“图像查看”-“射线图”，可以绘制射线图。图10-1射线图（2）单发射点实测曲线图点击“图像查看”-“实测曲线图（单点）”，可以绘制单发射点实测曲线图，可以利用工具栏中的按钮，实现不同发射点曲线的查看。图10-2单发射点实测曲线图（3）单发射