segh地震數(shù)據(jù)中浮點(diǎn)數(shù)的計(jì)算與分析

segh地震數(shù)據(jù)中浮點(diǎn)數(shù)的計(jì)算與分析

1浮點(diǎn)數(shù)的處理在不同的系統(tǒng)架構(gòu)下，浮點(diǎn)的數(shù)量在表示和存儲上存在較大差異。在SEGY地震資料中,數(shù)據(jù)可采用IEEE或者IBM格式浮點(diǎn)數(shù)存儲。特定計(jì)算機(jī)硬件平臺,只能處理與自己硬件平臺相同格式的浮點(diǎn)數(shù),若浮點(diǎn)數(shù)格式不屬于目標(biāo)處理機(jī)能夠處理的格式,則需進(jìn)行格式轉(zhuǎn)換。不同類型浮點(diǎn)數(shù)之間轉(zhuǎn)換,不僅涉及到階碼和尾數(shù)的表達(dá),還涉及到浮點(diǎn)數(shù)存儲的方式,以及位和位域的操作等問題。為此,本文設(shè)計(jì)了一個(gè)函數(shù)庫并提供API接口,用戶可以簡單直觀地實(shí)現(xiàn)兩種類型浮點(diǎn)數(shù)之間的轉(zhuǎn)換。232位置和沉沒點(diǎn)的形成結(jié)構(gòu)計(jì)算機(jī)中32位浮點(diǎn)數(shù)存儲在連續(xù)的4個(gè)字節(jié)中,不同體系架構(gòu)的計(jì)算機(jī)系統(tǒng),有不同的字節(jié)序。2.1處理器類型分析字節(jié)序是指多字節(jié)數(shù)據(jù)在計(jì)算機(jī)內(nèi)存中存儲或者網(wǎng)絡(luò)傳輸時(shí)各字節(jié)的存放順序,分為大字節(jié)序(BigEndian)、小字節(jié)序(LittleEndian)和混合序(Middle-endian)。每個(gè)數(shù)據(jù)都有所謂的“有效字節(jié)(Significantbyte)”,表示這個(gè)數(shù)據(jù)所占用的字節(jié)。例如一個(gè)32位整型數(shù),其有效位是4個(gè)字節(jié)。而對于0x0A0B0C0D來說,作為一個(gè)256進(jìn)制數(shù)(相對于10進(jìn)制數(shù)),其有效字節(jié)從高到低分別為0A、0B、0C及0D。大字節(jié)序(BigEndian)指其“最高有效字節(jié)(MostSignificantByte)”落在低地址上的存儲方式。例如地址a寫入0x0A0B0C0D之后,在內(nèi)存中的數(shù)據(jù):而小字節(jié)序(LittleEndian)正好相反,其“最低有效字節(jié)(LeastSignificantByte)”放在低地址上。例如:常用的CPU架構(gòu)如Intel、AMD的CPU使用的都是小字節(jié)序,而以前MacOS使用的是PowerPC大字節(jié)序(現(xiàn)在MacOS也使用Intel的CPU)。此外,還有一種很少見的混合序(middleendian),以2143的方式保存數(shù)據(jù)(相對于大字節(jié)序的1234及小字節(jié)序的4321),在此不再累述?？梢酝ㄟ^以下C函數(shù)來判斷處理器是大字節(jié)序還是小字節(jié)序:intcheckCPU(){intx=1;if(*(char*)&x==1)return1;elsereturn0;}若處理器是大字節(jié)序則返回0;若是小字節(jié)序,則返回1。下文敘述中,假設(shè)處理器類型都是小字節(jié)序。大字節(jié)序的處理器,可以經(jīng)過換位處理后進(jìn)行等效于小字節(jié)序的處理。2.2有效數(shù)字部分如表1浮點(diǎn)數(shù)是有理數(shù)中的特定子集,對于給定的實(shí)數(shù),在計(jì)算機(jī)中用近似值表示。具體來說,這個(gè)實(shí)數(shù)經(jīng)常用下式得到:V=(-1)S×M×A(E-B)(1)M=C+F(2)其中,V為實(shí)數(shù)值;S為符號,當(dāng)V為正數(shù)時(shí),取0,V為負(fù)數(shù)時(shí),取1;M為尾數(shù),代表實(shí)數(shù)值的有效數(shù)字部分,相當(dāng)于有效數(shù)字的整數(shù)部分(C)加上有小數(shù)部分(F);A為基數(shù);E為指數(shù);B為指數(shù)偏移;換句話說,這個(gè)實(shí)數(shù)由一個(gè)有效數(shù)字部分(即尾數(shù))乘以某個(gè)基數(shù)的整數(shù)次冪得到的,這種表示方法類似于基數(shù)為10的科學(xué)記數(shù)法,其中,根據(jù)計(jì)算機(jī)架構(gòu)(如IEEE和IBM架構(gòu))選用不同的常量A、B、C。2.3ieee浮點(diǎn)數(shù)的編碼32位IEEE單精度浮點(diǎn)數(shù)存儲在4個(gè)連續(xù)字節(jié)中,如圖3。32位IBM單精度浮點(diǎn)數(shù)存儲在4個(gè)連續(xù)字節(jié)中,如圖4。S代表符號,E(Exponent)代表指數(shù)部分,F(Fraction)代表小數(shù)部分。V=(-1)S×M×A(E-B)M=C+F32位IEEE浮點(diǎn)數(shù):A=2,B=127,C=1;32位IBM浮點(diǎn)數(shù):A=16,B=64,C=0。需要注意的是,無論是IEEE浮點(diǎn)數(shù)還是IBM浮點(diǎn)數(shù),都用了特定的32位二進(jìn)制數(shù)編碼分別表示0(包括+0和-0)、無窮大(包括正無窮大和負(fù)無窮大)和非數(shù)字。這些特定二進(jìn)制數(shù)編碼不符合浮點(diǎn)數(shù)基本計(jì)算公式,故在32位IEEE浮點(diǎn)數(shù)和IBM浮點(diǎn)數(shù)之間需要進(jìn)行對應(yīng)轉(zhuǎn)換,其對應(yīng)關(guān)系如表1。SEGY地震資料中,0經(jīng)常出現(xiàn),而無窮大和非數(shù)字較少出現(xiàn)。因此,本文提供的算法和代碼中,只考慮0在32位IEEE浮點(diǎn)數(shù)和IBM浮點(diǎn)數(shù)之間的轉(zhuǎn)換。2.43f個(gè)位為0時(shí)上文已經(jīng)提及,F區(qū)域用于存儲浮點(diǎn)數(shù)的小數(shù)部分。顯然,對于IEEE格式浮點(diǎn)數(shù),F的最小值為1/223(非0);對于IBM非0的浮點(diǎn)數(shù),F的最小值為1/16,而不是1/224。因?yàn)镕最多只能前3位為0,若第4位也為0,則違反了IBM格式的存儲規(guī)則。當(dāng)F個(gè)位為1時(shí),小數(shù)部分達(dá)到最大值。根據(jù)等比數(shù)列求和公式,可以得出:1/2+1/4+1/8+…+1/2n=1-1/2(n-1)以G代表小數(shù)區(qū)域的最大值,則:32位IEEE浮點(diǎn)數(shù):G=1-1/22232位IBM浮點(diǎn)數(shù):G=1-1/223由于浮點(diǎn)數(shù)有效數(shù)字部分M=C+F,可以算出各格式下浮點(diǎn)數(shù)的有效數(shù)字范圍:32位IEEE浮點(diǎn)數(shù):1<=M<=2-1/22232位IBM浮點(diǎn)數(shù):1/16<=M<=1-1/2232.53浮點(diǎn)數(shù)的選取以E代表浮點(diǎn)數(shù)的指數(shù)部分,實(shí)際值為E-B。其中B為指數(shù)偏移值,目的是實(shí)現(xiàn)負(fù)指數(shù)。32位IEEE浮點(diǎn)數(shù):0<=E<=254(255代表特定概念)32位IBM浮點(diǎn)數(shù):0<=E<=126(127代表特定概念)設(shè)T為浮點(diǎn)數(shù)絕對值的大小,則32位IEEE浮點(diǎn)數(shù):(1+1/223)×2-127<=T<=(2-1/222)×212732位IBM浮點(diǎn)數(shù):(1/16)×16-64<=T<=(1-1/223)×1662討論兩種格式浮點(diǎn)數(shù)的取值范圍是為了下文的溢出和舍入處理。負(fù)數(shù)則同理可得。332ieee單精度浮點(diǎn)數(shù)和vm單精度浮點(diǎn)數(shù)的轉(zhuǎn)化以V為實(shí)數(shù)值,32位IEEE浮點(diǎn)數(shù):V=(-1)S×(1+F)×2(E-127)(3)32位IBM浮點(diǎn)數(shù):V=(-1)S×F×16(E-64)(4)經(jīng)過簡單數(shù)學(xué)變換,可得IEEE單精度浮點(diǎn)數(shù)和IBM單精度浮點(diǎn)數(shù)的轉(zhuǎn)換規(guī)則。3.13生成領(lǐng)域小序或統(tǒng)收小序。請看浮點(diǎn)數(shù)字節(jié)序或成小節(jié)序,單位以小在進(jìn)行轉(zhuǎn)換之前,首先需要判斷計(jì)算機(jī)的類型(小字節(jié)序還是大字節(jié)序),并統(tǒng)一轉(zhuǎn)化為小字節(jié)序;其次設(shè)置掩碼與數(shù)據(jù)獲取相應(yīng)浮點(diǎn)數(shù)的S、E、F部分;然后再進(jìn)行上述計(jì)算;最后根據(jù)計(jì)算機(jī)類型轉(zhuǎn)化成相應(yīng)的字節(jié)序。3.1.1向度檢驗(yàn)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)上述公式中,數(shù)學(xué)符號「┐為向上取整。乘、除法用移位來實(shí)現(xiàn)。當(dāng)IEEE浮點(diǎn)數(shù)為特定概念時(shí),轉(zhuǎn)化為相應(yīng)IBM特定概念的二進(jìn)制編碼。EIEEE、FIEEE分別代表32位IEEE浮點(diǎn)數(shù)的指數(shù)、尾數(shù)部分;EIBM、FIBM代表32位IBM浮點(diǎn)數(shù)的指數(shù)、尾數(shù)部分。對于EIEEE,0<=EIEEE<=254,32<=(EIEEE+130)/4<=96,不會有上限溢出。但因?yàn)榇嬖谛?shù),需要FIBM進(jìn)行相應(yīng)的變換(除以相應(yīng)的數(shù))來獲得正確值,此時(shí)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,一般情況下不會有因求近似值而出現(xiàn)的誤差。但當(dāng)(EIEEE+130)/4余數(shù)為x時(shí),FIBM需要除以16(4-x)/4,即24-x。這時(shí)FIBM最后4-x(可能值為1,2,3)位將被舍去,造成誤差。3.1.2模型1:ieee浮點(diǎn)數(shù)對接枝型性能的影響每個(gè)浮點(diǎn)數(shù)用4個(gè)字節(jié)的unsignedlong類型存儲以方便轉(zhuǎn)換。unsignedlongIEEE2IBM_C(unsignedlongfloat_IEEE){if((float_IEEE<<1)==0){//若為IEEE浮點(diǎn)數(shù)//返回對應(yīng)特定IBM二進(jìn)制編碼returnfloat_IEEE;}//獲得IEEE浮點(diǎn)數(shù)的S,符號部分unsignedlongS_IEEE_32=float_IEEE&MASK_S;//獲得IEEE浮點(diǎn)數(shù)的E,指數(shù)部分unsignedlongE_IEEE_32=float_IEEE&MASK_IEEE32_E;//獲得IEEE浮點(diǎn)數(shù)的F,小數(shù)部分unsignedlongF_IEEE_32=float_IEEE&MASK_IEEE32_F;//獲得IBM浮點(diǎn)數(shù)的E,指數(shù)部分//E_IBM_32=(E_IEEE+130)/4//將E_IEEE放置在相應(yīng)的E_IBM的位置//獲取有效的指數(shù)部分的位數(shù)unsignedlongE_IBM_32=((E_IEEE_32+0x41000000)>>1)&MASK_IBM32_E;//獲得IBM浮點(diǎn)數(shù)的F,小數(shù)部分//將F_IEEE放置在相應(yīng)的F_IBM的位置//F_IBM=(F_IEEE+1)/2unsignedF_IBM_32=((F_IEEE_32<<1)+0x1000000)>>1;//但是當(dāng)(EIEEE+130)/4余數(shù)為x時(shí),//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浮點(diǎn)數(shù)16.5,-16.5,0,16.19(16.190001)轉(zhuǎn)化為IBM浮點(diǎn)數(shù),其對應(yīng)關(guān)系如表2。從表2可以看出,代碼將16.5,-16.5和0從IEEE格式浮點(diǎn)數(shù)轉(zhuǎn)換為IBM格式浮點(diǎn)數(shù)不會有任何誤差。但對于16.19,由于存在十進(jìn)制數(shù)和IEEE格式浮點(diǎn)數(shù)之間的轉(zhuǎn)換誤差,IEEE格式的浮點(diǎn)數(shù)只能用16.190001來近似表示16.19。同時(shí),IEEE格式浮點(diǎn)數(shù)16.190001最后的3位二進(jìn)制碼位均為1,其指數(shù)部分E=131,(131+130)/4余數(shù)為1,在轉(zhuǎn)換過程中最后3位被舍去,造成誤差,致使轉(zhuǎn)化后的IBM格式浮點(diǎn)數(shù)結(jié)果不為16.19。由此可見,上述代碼能夠成功將IEEE格式浮點(diǎn)數(shù)轉(zhuǎn)化為IBM格式浮點(diǎn)數(shù)。3.23兩位ibm浮點(diǎn)的方法和代碼轉(zhuǎn)換為ee浮點(diǎn)的數(shù)量3.2.1fie型浮點(diǎn)數(shù)EIEEE=4EIBM-130(7)FIEEE=2FIBM-1(8)(其中,當(dāng)IBM浮點(diǎn)數(shù)為特定概念時(shí),轉(zhuǎn)化為相應(yīng)IEEE特定概念的二進(jìn)制編碼)對于FIBM,1/16<=FIBM<=1-1/223,1/8-1<=(2FIBM-1)<=1-1/222顯然,當(dāng)FIBM<1/2,FIEEE將溢出,此時(shí)需要進(jìn)行調(diào)整。由于FIBM>=1/16,設(shè)第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浮點(diǎn)數(shù)以16為基底,浮點(diǎn)數(shù)的大小范圍(1/16)×16-64<=T<=(1-1/223)×1662大于IEEE浮點(diǎn)數(shù)。因此,只有當(dāng)0<=4EIBM-130<=254,即32<=EIBM<=96時(shí),才能正確地將IBM浮點(diǎn)數(shù)轉(zhuǎn)為IEEE浮點(diǎn)數(shù)。當(dāng)EIBM<32時(shí),轉(zhuǎn)換為IEEE格式浮點(diǎn)數(shù)結(jié)果為0;當(dāng)EIBM>96時(shí),轉(zhuǎn)換為IBM格式浮點(diǎn)數(shù)結(jié)果為無窮大。當(dāng)FIBM前3位為0時(shí),EIBM=97,其IBM浮點(diǎn)數(shù)一樣可以正確地轉(zhuǎn)化為IEEE浮點(diǎn)數(shù)。3.2.2mask現(xiàn)代浮點(diǎn)數(shù)的編碼其中每個(gè)浮點(diǎn)數(shù)用4個(gè)字節(jié)的unsignedlong類型存儲,方便轉(zhuǎn)換unsignedlongIBM2IEEE_C(unsignedlongfloat_IBM){if((float_IBM<<1)==0){//若為IBM浮點(diǎn)數(shù)//返回對應(yīng)特定IEEE二進(jìn)制編碼returnfloat_IBM;}//獲得IBM浮點(diǎn)數(shù)的S,符號部分unsignedlongS_IBM_32=float_IBM&MASK_S;//獲得IBM浮點(diǎn)數(shù)的E,指數(shù)部分unsignedlongE_IBM_32=float_IBM&MASK_IBM32_E;//獲得IBM浮點(diǎn)數(shù)的F,小數(shù)部分unsignedlongF_IBM_32=float_IBM&MASK_IBM32_F;//獲得IBM浮點(diǎn)數(shù)的F,小數(shù)部分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部分的相應(yīng)位置F_IEEE_32=(F_IEEE_32-0x01000000)>>1;//獲得IBM浮點(diǎn)數(shù)的E,指數(shù)部分//開始計(jì)算//放置在正確位置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浮點(diǎn)數(shù)16.5,-16.5,0,16.19(16.18998718)轉(zhuǎn)化為IEEE浮點(diǎn)數(shù)格式,其對應(yīng)關(guān)系如表3。從表3可以看出,代碼將16.5,-16.5和0從IBM格式浮點(diǎn)數(shù)轉(zhuǎn)換為IEEE格式浮點(diǎn)數(shù)不會有任何誤差。但是對于16.19,由于存在十進(jìn)制數(shù)和IBM格式浮點(diǎn)數(shù)之間的轉(zhuǎn)換誤差,IBM格式的浮點(diǎn)數(shù)只能用16.18998718來近似表示16.19,與直接將16.19用IEEE格式浮點(diǎn)數(shù)存儲存在差別。由于IBM格式浮點(diǎn)數(shù)轉(zhuǎn)換后得到的IEEE格式浮點(diǎn)數(shù)編碼(0x41818518)和原來的IEEE浮點(diǎn)數(shù)編碼(0x4181851F)的差別僅限于最后3位二進(jìn)制數(shù),其中一個(gè)是“0”,另一個(gè)為“1”,說明差別是由轉(zhuǎn)換過程中的舍去造成,而不是程序錯(cuò)誤。由此可見,上述代碼能夠成功地將IBM格式浮點(diǎn)數(shù)轉(zhuǎn)化為IEEE格式浮點(diǎn)數(shù)。值得注意的是,IBM格式的浮點(diǎn)數(shù)由于以16為基數(shù),能表示較大的數(shù)據(jù)范圍,其尾數(shù)雖然是24位二進(jìn)制數(shù),但是只能保證21位的有效位,精度上較IEEE格式的浮點(diǎn)數(shù)要低。當(dāng)把IEEE浮點(diǎn)轉(zhuǎn)換成IBM浮點(diǎn)數(shù)時(shí),會出現(xiàn)精度降低的現(xiàn)象,而在把IBM浮點(diǎn)數(shù)轉(zhuǎn)換為IEEE浮點(diǎn)數(shù)時(shí),會出現(xiàn)溢出現(xiàn)象。因此,其轉(zhuǎn)換只是最大限度的近似,不完全可逆。4保證單次存儲float信息在SEGY格式的地震資料中,數(shù)據(jù)以道為單位存放。地震道數(shù)據(jù)由道頭段和數(shù)據(jù)段組成,道頭段中有浮點(diǎn)數(shù)格式等參數(shù)信息,數(shù)據(jù)段由包含n個(gè)樣值點(diǎn)的數(shù)據(jù)組成。在實(shí)際處理中,其浮點(diǎn)數(shù)用一維數(shù)組的數(shù)據(jù)結(jié)構(gòu)進(jìn)行存儲,形成一個(gè)輸入的浮點(diǎn)數(shù)數(shù)組buffer。這時(shí),可以利用下面函數(shù),將整個(gè)buffer浮點(diǎn)數(shù)數(shù)組轉(zhuǎn)換成另一種浮點(diǎn)數(shù)。源代碼如下:(其中TYPE為轉(zhuǎn)換方向參數(shù),IEEE轉(zhuǎn)IBM則TYPE為0;IBM轉(zhuǎn)IEEE則TYPE為1。n為地震數(shù)據(jù)的樣值數(shù),buffer為輸入輸出緩沖區(qū)。)intfloatConversion(float