談談體數據2：體數據格式轉換工具vodaconv的使用及原理

談談體數據2：體數據格式轉換工具vodaconv的使用及原理

文/Sobereva @北京科音
First release: 2012-Feb-14  Last update: 2012-Feb-22

前言：在《談談體數據1：介紹體數據》（http://www.shanxitv.org/127）一文中寡人介紹了體數據的相關概念。raw是廣為通用的體數據儲存格式，可惜計算化學領域的可視化軟件極少有支持raw格式的；而通用的渲染體數據的程序也沒有能支持計算化學領域通用的Gaussian型cube文件(.cub)的。為了能夠讓體數據互用，架起領域間的橋梁，寡人開發了vodaconv程序用于格式轉換，將在本文介紹。在《談談體數據》系列的接下來的文章中將會用到這個程序對raw和cube型體數據進行相互轉換。同時，本文也將仔細介紹此程序的原理，將有助于讀者了解與體數據相關的一些細節問題。

1 vodaconv程序介紹

體數據轉換工具vodaconv支持的輸入格式包括raw、Gaussian型cube文件(.cub)、OpenDX的格點數據格式(.dx)、Dmol3的.grd格式。支持的輸出格式包括raw和cube。在未來的版本中將會支持更多格式。此程序只支持立方格子的格點文件，即平移向量平行于笛卡爾坐標軸。

此程序可以在此下載：/usr/uploads/file/20150611/20150611014131_45262.rar。

程序會根據擴展名判斷格式。如果輸入的是.raw文件，需要由用戶自行輸入各個方向的格點數和格點間隔、起始位置以及數據類型。如果輸入不對的話將導致此程序在讀取.raw文件中自動退出。可以事先自行算一下輸入的格點數和數據類型與當前文件是否相符，比如含有256*256*128個格點的16bit的raw文件，文件大小應當為256*256*128*16/8=16777216字節，如果文件大小與之不符，說明數據類型或者格點數與此文件不符。

程序目前支持的raw數據類型包括8bit整數、16bit整數和單精度浮點數三種。讀入、輸出raw文件時整數一律按照無符號整數形式處理，即8bit整數數據范圍為0~255，16bit整數數據范圍為0~65535。vodaconv.exe讀寫的raw文件是按照Big Endian順序，而vodaconv_LE.exe是按照Little Endian順序（含義將在后文解釋）。

在讀入格點數據后會提示是否對數據進行映射，需要先輸入原先的數據范圍上下限，超過上限的數據將被設為上限值，超過下限的數據將被設為下限值。然后需要輸入映射后的數據范圍上下限。原先的數據將按照上下限范圍的比例映射到新的數值范圍內。這個功能主要用于調節數值范圍，以免以8/16bit整數方式輸出raw文件時數據越界。利用這個功能也可以將16bit整數的raw格式轉化為8bit整數的raw格式，即downsample過程。一個應用例子是，比如想把電子定域化函數(ELF)的cube文件轉化成16bit整數raw文件，由于已知ELF函數值域是[0,1]，原先數據范圍就可以輸入0,1，而映射后的數據范圍就可以輸入0,65535，這樣ELF函數值域區間就均勻展開到了無符號16bit整數的記錄范圍里了。

之后會提示用戶輸入要輸出的文件名，也是根據后綴名判斷輸出格式。在輸出cube文件時，若之前載入的是.raw/.dx/.grd格式，由于其中沒有定義原子信息，為了避免一些程序在讀取cube文件時必須要求至少有一個原子造成的不兼容，此程序將自動在cube文件內添加一個位于原點的氫。

2 vodaconv程序的編譯方法

程序壓縮包內已包含win32平臺下由Intel visual fortran 12.0預編譯好的可執行文件。同時附帶了源代碼文件。

若想自行編譯源代碼，在Linux下執行ifort vodaconv.f90 -o vodaconv即可。不要妄圖用gfortran編譯，因為此程序用了integer*1的數據類型，gfortran不支持。

如果想在Windows下編譯，用Compaq visual fortran或Intel visual fortran都可以。首先將vodaconv.f90拖進新建的項目里，然后設置stack size為較大值，比如256MB（否則由于默認棧內存過小，在載入較大的格點文件時會導致棧內存溢出），之后照常編譯即可。在CVF中設置stack size的方法是選擇"Project"-"Settings..."-"Link"，在"Category"里選擇"Output"，在"Reserve"里輸入"0x10000000"（即256MB，這也是CVF支持的上限）。在IVF中設置stack size的方法是選擇“項目”-“XX屬性”，選"Linker"，將"System"里的"Stack Reserve Size"設為"268435456"（即256MB。當然也可以設得更大使程序能載入巨型格點數據）。

如果你希望你所編譯出的本程序使用Little Endian順序讀寫二進制文件，那么就按上面的步驟編譯就行了，編譯器默認就是用Little Endian。如果希望使用Big Endian順序，那么在Linux下編譯時應該加上-convert big_endian選項。如果是在IVF中，應當在項目的屬性的Fortran-Compatibility中將Unformatted file conversion設成Big Endian。如果是在CVF中，進入"Project"-"Settings-"Fortran"，選Compatibility分類，將Unformatted file conversion設成Big Endian。

3 關于含符號整數和無符號整數

這里介紹一下含符號和無符號整數與raw文件的關系。這里討論的都是Fortran語言涉及的整數。

計算機程序中用到的整數數據一般是雙字節整數（16bit整數。也叫短整數，short integer）、四字節整數（32bit整數。也叫長整數，long integer），在現代的編譯器中（比如版本較新的ifort）還支持八字節整數。單字節整數，即8bit整數一般很少用到，因為它的數值范圍太窄，一般沒什么用。而且它并不是Fortran標準中定義的數據類型，但CVF、ifort等編譯器是支持它的（用integer*1定義。在gfortran中不支持）。然而8bit整數卻被普遍被raw文件用于記錄低精度數據，因此我們不能忽視它。

無符號整數意思是整數值只為正，因此8bit無符號整數的范圍是0~2^8-1，即0~255；16bit無符號整數是0~2^16-1，即0~65535。含符號整數則將一半記錄空間用于記錄負值，8bit和16bit含符號整數的數值范圍分別是-128~127和-32768~32767。

在Fortran標準中只定義了含符號整數，然而raw格式中的整數數據都是使用無符號整數來記錄的，所以必須搞清楚含符號和無符號整數間的轉換關系，才能用Fortran程序正確載入整數型raw格式數據。

在計算機的內存或二進制文件中，四個十六進制數由小到大對應的16bit整數數值如下所示，括號中第一個對應的是含符號型整數值，第二個是無符號型整數值：
0x0000(0,0)->0x0001(1,1)->...->0x7FFF(32767,32767)->0x8000(-32768,32768)->0x8001(-32767,32769)->...->0xFFFF(-1,65535)。
即曰：對于0~32767的部分，含符號和無符號的整數是以相同方式記錄的。對于含符號整型，如果大于了32767就會成為負值，此時這個負值加上65536就可以變為相應的無符號整型的數值。這個規則對于其它長度的整型數據也是一樣。

由于Fortran用的都是含符號整型，這也就是說，raw文件記錄的8bit/16bit無符號整型數據會被當作是含符號整型數據被讀入，因此需要在讀入后把所有負值數據都加上256（對于8bit）和65536（對于16bit），這樣才是正確的數值。

在通過write語句輸出8bit和16bit整數型raw文件前，需要先分別利用int1()和int2()函數將內存中的數據轉換為8bit和16bit含符號整型。對于8bit整型，大于127的數據將被自動減去256成為負數使之能被含符號整型變量表示；對于16bit則會自動減去65536。這個過程和載入時需要做的變換是互逆的。

4 二進制文件的字節記錄順序問題

不同的平臺、不同編譯器設定下產生的不同程序所利用、生成的二進制文件中的多字節數據有不同的字節記錄順序。8bit整數就一個字節所以無需考慮此問題，而16bit整數有兩個字節，就會牽扯到這個問題。簡要來說，例如32444這個16bit整數對應的十六進制數是0x7EBC，這個數據在二進制文件中如果低字節記錄的是0xBC，高字節是0x7E，那么就稱為Little Endian（小端）順序，用Ultraedit等支持十六進制的編輯器打開此文件就會看到相應位置顯示的是BC 7E。若低、高字節分別記錄的是0x7E和0xBC，就稱為Big Endian（大端）順序。Windows、Linux等系統下的程序通常用的是Little Endian，而據說Mac OS用的是Big Endian。

如果你發現用vodaconv.exe（用的是Big Endian）轉換的16bit整數型raw文件在載入其它一些體數據可視化程序后數值不正確，那么應該嘗試用vodaconv_LE.exe（用的是Little Endian）重新生成raw文件。而一些可視化程序自身同時支持多種字節記錄順序，可在載入文件時選擇。

5 vodaconv程序的原理

此程序源代碼很容易讀懂，沒必要有太多注釋。名為defvar的module用于保存格點設定。orgx,orgy,orgz是原點三個分量，dx,dy,dz是三個方向格點平移距離，nx,ny,nz是三個方向的點數。數組a保存原子信息，包括坐標，電荷，原子序數。cubmati1、cubmati2和cubmatr4分別是8bit含符號整型、16bit含符號整型和單精度浮點型動態分配空間的數組，用于儲存格點數據。

下面介紹一下程序流程

1 載入體數據文件：首先根據輸入文件的擴展名判定輸入文件的格式。由于raw文件不含格點設定信息，而且數據類型未知，故需要由用戶輸入。平移矢量和原點坐標可輸入可不輸入，按ENTER可以使用默認值，這會將坐標原點設定到體數據的中心位置。.dx、.grd和.cub文件格式大同小異，載入.cub文件的代碼在《Gaussian型cube文件簡介及讀、寫方法和簡單應用》（http://www.shanxitv.org/125）一文中已經仔細介紹過了。這三種格式都是以浮點數表示數據，由于數據精度都不很高，為節約內存，用單精度浮點數記錄就足夠了。載入的數據暫存在cubmatr4數組中。

對于raw格式，若是單精度浮點型就直接載入到cubmatr4數組中。如果是8bit和16bit無符號整型，就分別載入到cubmati1和cubmati2中，然后將其數值轉移到cubmatr4中，然后根據第三節所述對cubmatr4的負值數據部分進行運算以還原成原本的數值。另外，這也就是說在vodaconv程序內部對任何數據類型一律以單精度浮點方式保存，因為單精度浮點的儲存精度和數值范圍都遠遠高于8/16bit整型。

載入完畢后會顯示體數據的格點設定信息以及最大和最小值。

2 對載入的數據的數值范圍進行變換：這是給用戶提供一個便利，可以避免數據越界，也可以將某一數值區域進行“放大”以使感興趣的數值范圍容易展現。在第一節也已經討論了。如果想對數據做更豐富、靈活的變換操作，建議用寡人編寫的Multiwfn程序2.3及以后版本（http://www.shanxitv.org/multiwfn），它的主功能13專門用來操縱和提取cube格式的格點數據。

3 輸出體數據文件：根據后綴名判斷要輸出的文件類型，將cubmatr4中的數據以相應方式輸出即可。輸出cube文件的代碼在《Gaussian型cube文件簡介及讀、寫方法和簡單應用》里也已經詳細討論了。如本文第一節所述，可能會在cube文件中加入一個氫原子信息。