segh地震数据中浮点数的计算与分析

segh地震数据中浮点数的计算与分析

1浮点数的处理在不同的系统架构下，浮点的数量在表示和存储上存在较大差异。在SEGY地震资料中,数据可采用IEEE或者IBM格式浮点数存储。特定计算机硬件平台,只能处理与自己硬件平台相同格式的浮点数,若浮点数格式不属于目标处理机能够处理的格式,则需进行格式转换。不同类型浮点数之间转换,不仅涉及到阶码和尾数的表达,还涉及到浮点数存储的方式,以及位和位域的操作等问题。为此,本文设计了一个函数库并提供API接口,用户可以简单直观地实现两种类型浮点数之间的转换。232位置和沉没点的形成结构计算机中32位浮点数存储在连续的4个字节中,不同体系架构的计算机系统,有不同的字节序。2.1处理器类型分析字节序是指多字节数据在计算机内存中存储或者网络传输时各字节的存放顺序,分为大字节序(BigEndian)、小字节序(LittleEndian)和混合序(Middle-endian)。每个数据都有所谓的“有效字节(Significantbyte)”,表示这个数据所占用的字节。例如一个32位整型数,其有效位是4个字节。而对于0x0A0B0C0D来说,作为一个256进制数(相对于10进制数),其有效字节从高到低分别为0A、0B、0C及0D。大字节序(BigEndian)指其“最高有效字节(MostSignificantByte)”落在低地址上的存储方式。例如地址a写入0x0A0B0C0D之后,在内存中的数据:而小字节序(LittleEndian)正好相反,其“最低有效字节(LeastSignificantByte)”放在低地址上。例如:常用的CPU架构如Intel、AMD的CPU使用的都是小字节序,而以前MacOS使用的是PowerPC大字节序(现在MacOS也使用Intel的CPU)。此外,还有一种很少见的混合序(middleendian),以2143的方式保存数据(相对于大字节序的1234及小字节序的4321),在此不再累述。可以通过以下C函数来判断处理器是大字节序还是小字节序:intcheckCPU(){intx=1;if(*(char*)&x==1)return1;elsereturn0;}若处理器是大字节序则返回0;若是小字节序,则返回1。下文叙述中,假设处理器类型都是小字节序。大字节序的处理器,可以经过换位处理后进行等效于小字节序的处理。2.2有效数字部分如表1浮点数是有理数中的特定子集,对于给定的实数,在计算机中用近似值表示。具体来说,这个实数经常用下式得到:V=(-1)S×M×A(E-B)(1)M=C+F(2)其中,V为实数值;S为符号,当V为正数时,取0,V为负数时,取1;M为尾数,代表实数值的有效数字部分,相当于有效数字的整数部分(C)加上有小数部分(F);A为基数;E为指数;B为指数偏移;换句话说,这个实数由一个有效数字部分(即尾数)乘以某个基数的整数次幂得到的,这种表示方法类似于基数为10的科学记数法,其中,根据计算机架构(如IEEE和IBM架构)选用不同的常量A、B、C。2.3ieee浮点数的编码32位IEEE单精度浮点数存储在4个连续字节中,如图3。32位IBM单精度浮点数存储在4个连续字节中,如图4。S代表符号,E(Exponent)代表指数部分,F(Fraction)代表小数部分。V=(-1)S×M×A(E-B)M=C+F32位IEEE浮点数:A=2,B=127,C=1;32位IBM浮点数:A=16,B=64,C=0。需要注意的是,无论是IEEE浮点数还是IBM浮点数,都用了特定的32位二进制数编码分别表示0(包括+0和-0)、无穷大(包括正无穷大和负无穷大)和非数字。这些特定二进制数编码不符合浮点数基本计算公式,故在32位IEEE浮点数和IBM浮点数之间需要进行对应转换,其对应关系如表1。SEGY地震资料中,0经常出现,而无穷大和非数字较少出现。因此,本文提供的算法和代码中,只考虑0在32位IEEE浮点数和IBM浮点数之间的转换。2.43f个位为0时上文已经提及,F区域用于存储浮点数的小数部分。显然,对于IEEE格式浮点数,F的最小值为1/223(非0);对于IBM非0的浮点数,F的最小值为1/16,而不是1/224。因为F最多只能前3位为0,若第4位也为0,则违反了IBM格式的存储规则。当F个位为1时,小数部分达到最大值。根据等比数列求和公式,可以得出:1/2+1/4+1/8+…+1/2n=1-1/2(n-1)以G代表小数区域的最大值,则:32位IEEE浮点数:G=1-1/22232位IBM浮点数:G=1-1/223由于浮点数有效数字部分M=C+F,可以算出各格式下浮点数的有效数字范围:32位IEEE浮点数:1<=M<=2-1/22232位IBM浮点数:1/16<=M<=1-1/2232.53浮点数的选取以E代表浮点数的指数部分,实际值为E-B。其中B为指数偏移值,目的是实现负指数。32位IEEE浮点数:0<=E<=254(255代表特定概念)32位IBM浮点数:0<=E<=126(127代表特定概念)设T为浮点数绝对值的大小,则32位IEEE浮点数:(1+1/223)×2-127<=T<=(2-1/222)×212732位IBM浮点数:(1/16)×16-64<=T<=(1-1/223)×1662讨论两种格式浮点数的取值范围是为了下文的溢出和舍入处理。负数则同理可得。332ieee单精度浮点数和vm单精度浮点数的转化以V为实数值,32位IEEE浮点数:V=(-1)S×(1+F)×2(E-127)(3)32位IBM浮点数:V=(-1)S×F×16(E-64)(4)经过简单数学变换,可得IEEE单精度浮点数和IBM单精度浮点数的转换规则。3.13生成领域小序或统收小序。请看浮点数字节序或成小节序,单位以小在进行转换之前,首先需要判断计算机的类型(小字节序还是大字节序),并统一转化为小字节序;其次设置掩码与数据获取相应浮点数的S、E、F部分;然后再进行上述计算;最后根据计算机类型转化成相应的字节序。3.1.1向度检验EIBM=|(EIEEE+130)/4|(5)EΙBΜ=|(EΙEEE+130)/4|(5)FIBM=(FIEEE+1)2×16⌈(EIEEE+130)/4⌉−(EIEEE+130)/4(6)FΙBΜ=(FΙEEE+1)2×16⌈(EΙEEE+130)/4⌉-(EΙEEE+130)/4(6)上述公式中,数学符号「┐为向上取整。乘、除法用移位来实现。当IEEE浮点数为特定概念时,转化为相应IBM特定概念的二进制编码。EIEEE、FIEEE分别代表32位IEEE浮点数的指数、尾数部分;EIBM、FIBM代表32位IBM浮点数的指数、尾数部分。对于EIEEE,0<=EIEEE<=254,32<=(EIEEE+130)/4<=96,不会有上限溢出。但因为存在小数,需要FIBM进行相应的变换(除以相应的数)来获得正确值,此时EIEEE向上取整。FIBM×16×EIBM=C即EIBM每增加x,FIBM要除以16x。其中x的可能取值为0.25,0.5,0.75。对于FIBM,(1+1/223)<=(FIEEE+1)<=2-1/222,1/2+1/224<=(FIEEE+1)/2<=1-1/223,一般情况下不会有因求近似值而出现的误差。但当(EIEEE+130)/4余数为x时,FIBM需要除以16(4-x)/4,即24-x。这时FIBM最后4-x(可能值为1,2,3)位将被舍去,造成误差。3.1.2模型1:ieee浮点数对接枝型性能的影响每个浮点数用4个字节的unsignedlong类型存储以方便转换。unsignedlongIEEE2IBM_C(unsignedlongfloat_IEEE){if((float_IEEE<<1)==0){//若为IEEE浮点数//返回对应特定IBM二进制编码returnfloat_IEEE;}//获得IEEE浮点数的S,符号部分unsignedlongS_IEEE_32=float_IEEE&MASK_S;//获得IEEE浮点数的E,指数部分unsignedlongE_IEEE_32=float_IEEE&MASK_IEEE32_E;//获得IEEE浮点数的F,小数部分unsignedlongF_IEEE_32=float_IEEE&MASK_IEEE32_F;//获得IBM浮点数的E,指数部分//E_IBM_32=(E_IEEE+130)/4//将E_IEEE放置在相应的E_IBM的位置//获取有效的指数部分的位数unsignedlongE_IBM_32=((E_IEEE_32+0x41000000)>>1)&MASK_IBM32_E;//获得IBM浮点数的F,小数部分//将F_IEEE放置在相应的F_IBM的位置//F_IBM=(F_IEEE+1)/2unsignedF_IBM_32=((F_IEEE_32<<1)+0x1000000)>>1;//但是当(EIEEE+130)/4余数为x时,//FIBM要除以xunsignedlongreminder=(E_IEEE_32+0x41000000)>>23;if((reminder&3)!=0){E_IBM_32=E_IBM_32+0x01000000;F_IBM_32=F_IBM_32>>(4-(reminder&3));}returnS_IEEE_32|E_IBM_32|F_IBM_32;}32位IEEE浮点数16.5,-16.5,0,16.19(16.190001)转化为IBM浮点数,其对应关系如表2。从表2可以看出,代码将16.5,-16.5和0从IEEE格式浮点数转换为IBM格式浮点数不会有任何误差。但对于16.19,由于存在十进制数和IEEE格式浮点数之间的转换误差,IEEE格式的浮点数只能用16.190001来近似表示16.19。同时,IEEE格式浮点数16.190001最后的3位二进制码位均为1,其指数部分E=131,(131+130)/4余数为1,在转换过程中最后3位被舍去,造成误差,致使转化后的IBM格式浮点数结果不为16.19。由此可见,上述代码能够成功将IEEE格式浮点数转化为IBM格式浮点数。3.23两位ibm浮点的方法和代码转换为ee浮点的数量3.2.1fie型浮点数EIEEE=4EIBM-130(7)FIEEE=2FIBM-1(8)(其中,当IBM浮点数为特定概念时,转化为相应IEEE特定概念的二进制编码)对于FIBM,1/16<=FIBM<=1-1/223,1/8-1<=(2FIBM-1)<=1-1/222显然,当FIBM<1/2,FIEEE将溢出,此时需要进行调整。由于FIBM>=1/16,设第n位(可能取值为1,2,3,4)为FIBM中为1的最高位,则:EIEEE=(4EIBM-130)-(n-1)(9)FIEEE=2nFIBM-1(10)对于EIBM,由于0<=EIBM<=126,-130<=4EIBM-130<=374。IBM浮点数以16为基底,浮点数的大小范围(1/16)×16-64<=T<=(1-1/223)×1662大于IEEE浮点数。因此,只有当0<=4EIBM-130<=254,即32<=EIBM<=96时,才能正确地将IBM浮点数转为IEEE浮点数。当EIBM<32时,转换为IEEE格式浮点数结果为0;当EIBM>96时,转换为IBM格式浮点数结果为无穷大。当FIBM前3位为0时,EIBM=97,其IBM浮点数一样可以正确地转化为IEEE浮点数。3.2.2mask现代浮点数的编码其中每个浮点数用4个字节的unsignedlong类型存储,方便转换unsignedlongIBM2IEEE_C(unsignedlongfloat_IBM){if((float_IBM<<1)==0){//若为IBM浮点数//返回对应特定IEEE二进制编码returnfloat_IBM;}//获得IBM浮点数的S,符号部分unsignedlongS_IBM_32=float_IBM&MASK_S;//获得IBM浮点数的E,指数部分unsignedlongE_IBM_32=float_IBM&MASK_IBM32_E;//获得IBM浮点数的F,小数部分unsignedlongF_IBM_32=float_IBM&MASK_IBM32_F;//获得IBM浮点数的F,小数部分unsignedlongradix=0;unsignedlongF_IEEE_32=F_IBM_32;while(radix<=3&&F_IEEE_32<0x01000000){radix++;F_IEEE_32=F_IEEE_32<<1;}//放置回IEEE类型F部分的相应位置F_IEEE_32=(F_IEEE_32-0x01000000)>>1;//获得IBM浮点数的E,指数部分//开始计算//放置在正确位置unsignedlongE_IEEE_32=(((E_IBM_32>>22)-130)-(radix-1))<<23;//是否溢出if(E_IEEE_32>0x7F800000){returnS_IBM_32|0x7F800000;}if(E_IEEE_32<0x10000000){returnS_IBM_32;}elsereturnS_IBM_32|E_IEEE_32|F_IEEE_32;}32位IBM浮点数16.5,-16.5,0,16.19(16.18998718)转化为IEEE浮点数格式,其对应关系如表3。从表3可以看出,代码将16.5,-16.5和0从IBM格式浮点数转换为IEEE格式浮点数不会有任何误差。但是对于16.19,由于存在十进制数和IBM格式浮点数之间的转换误差,IBM格式的浮点数只能用16.18998718来近似表示16.19,与直接将16.19用IEEE格式浮点数存储存在差别。由于IBM格式浮点数转换后得到的IEEE格式浮点数编码(0x41818518)和原来的IEEE浮点数编码(0x4181851F)的差别仅限于最后3位二进制数,其中一个是“0”,另一个为“1”,说明差别是由转换过程中的舍去造成,而不是程序错误。由此可见,上述代码能够成功地将IBM格式浮点数转化为IEEE格式浮点数。值得注意的是,IBM格式的浮点数由于以16为基数,能表示较大的数据范围,其尾数虽然是24位二进制数,但是只能保证21位的有效位,精度上较IEEE格式的浮点数要低。当把IEEE浮点转换成IBM浮点数时,会出现精度降低的现象,而在把IBM浮点数转换为IEEE浮点数时,会出现溢出现象。因此,其转换只是最大限度的近似,不完全可逆。4保证单次存储float信息在SEGY格式的地震资料中,数据以道为单位存放。地震道数据由道头段和数据段组成,道头段中有浮点数格式等参数信息,数据段由包含n个样值点的数据组成。在实际处理中,其浮点数用一维数组的数据结构进行存储,形成一个输入的浮点数数组buffer。这时,可以利用下面函数,将整个buffer浮点数数组转换成另一种浮点数。源代码如下:(其中TYPE为转换方向参数,IEEE转IBM则TYPE为0;IBM转IEEE则TYPE为1。n为地震数据的样值数,buffer为输入输出缓冲区。)intfloatConversion(float