OSGB格式三維模型投影坐標轉換方法

程曉光，王婉秋，陸泉源

(飛燕航空遙感技術有限公司，南京 210001)

目前，基于傾斜攝影測量技術生產的三維模型已成為主要的測繪產品之一，其生產步驟包括密集點云提取、網格面模型構建、紋理映射、成果編輯、成果導出等[1]，較為復雜。ContextCapture[2]是目前傾斜攝影測量領域主要的三維建模軟件之一[1, 3]，而OpenSceneGraph Binary(OSGB)是ContextCapture導出三維模型時最常用的格式之一，由開源三維圖形工具包OpenSceneGraph[4](OSG)推出，屬于二進制文件格式，帶有嵌入式鏈接紋理數據，允許以細節層次(Level Of Detail，LOD)技術進行多級金字塔顯示[4]。

三維模型使用的空間參考系統以投影坐標系居多。根據國家標準[5-8]，在基本比例尺地形圖和正射影像地圖中，大地基準采用2000國家大地坐標系，高程基準主要采用1985國家高程基準，除1∶1 000 000外的其它比例尺的投影方法均為3°或6°分帶的高斯-克呂格投影。此外，各大城市和地區還建立了地方獨立坐標系，目的是將高程歸化與投影變形的影響控制在微小范圍內，保證城市范圍內長度變形最小，以滿足當地規劃建設需要[9]。大多數地方獨立坐標系是在1954年北京坐標系或1980年西安坐標系的基礎上利用高斯-克呂格投影進行改化后建立的。常見改化包括將中央子午線移至測區中央，將投影面改為測區平均高程面、對原點縱橫坐標進行偏移等[9]。在三維模型交付時，除了國家標準投影坐標系外，還經常需要以地方獨立坐標系交付。

使用三維建模軟件在不同投影坐標系下各生產一遍三維模型會造成巨大的工作量，且模型間可能存在偏差。如果ContextCapture在輸出模型時選擇的投影坐標系與空三使用的坐標系基于不同橢球體，則輸出模型的坐標易有較大誤差。有商業軟件可對三維模型進行投影坐標轉換，但轉換后的模型經常存在紋理丟失、顯示閃爍等問題。

針對ContextCapture導出的OSGB格式的三維模型的投影坐標轉換問題，文中分析OSGB格式的三維模型的目錄和文件結構，提出一種簡單高效的投影坐標轉換方法，先使用OSG獲得待轉換坐標，再進行坐標轉換，最后利用文件坐標數據替換法保存三維模型，并給出示例代碼。該方法擴展性強，轉換后的模型無紋理丟失、無變形、無閃爍。

1 OSGB格式模型的目錄和文件結構

在ContextCapture的生產項目定義中設置輸出格式為OSGB，劃分瓦片(Tile)，使用LOD技術，則輸出的三維模型存儲文件夾的結構見圖1。Data文件夾下為各瓦片的子文件夾，每個子文件夾存儲該瓦片不同金字塔級別的OSGB文件，其中一個OSGB文件和子文件夾同名。但是該模式并不便于模型瀏覽。為此，一種常見的做法是在三維模型存儲文件夾下，額外放置S3C文件以對構成場景的瓦片文件進行索引，方便瀏覽顯示。

圖1 OSGB格式的三維模型的目錄結構

為降低內存和磁盤使用，OSGB文件常采用4字節長度(單精度)而非8字節長度(雙精度)的浮點數來存儲頂點坐標。對中國境內的高斯-克呂格投影來說，如果坐標以m為單位，則整數部分常會有6位或7位，但單精度浮點數只能保證7位十進制有效數字[10]。如果在OSGB文件中存儲頂點原始投影坐標，則容易造成坐標小數部分有效位數不夠。為此，在采用投影坐標系的情況下，OSGB文件一般對頂點坐標做偏移，以保證小數點后有足夠的有效數字。坐標偏移量存儲在metadata.xml的SRSOrigin標簽中。在metadata.xml中，投影坐標系信息存儲在SRS標簽，支持歐洲石油調查組織[11](European Petroleum Survey Group，EPSG)定義的坐標系編碼或眾知文本(Well-Known Text，WKT)。

單個OSGB文件存儲金字塔某級在特定空間范圍內的模型，其基于OSG體系的典型組織結構見圖2。下文如無特殊說明，則OSG類的默認命名空間是osg。位于金字塔最底層的模型的根節點為一個Geode類對象；利用四叉樹或八叉樹等[12]技術做模型分割的根節點為一個Group類對象；除最底層外，未做分割的模型的根節點為一個PagedLOD類對象。這其中涉及到的主要節點類之間的繼承關系見圖3。

圖2 典型OSGB文件的組織結構

圖3 主要節點類之間的繼承關系(箭頭起點為父類，終點為子類)

2 投影坐標轉換方法

對OSGB格式的三維模型進行投影坐標轉換，可以先對模型頂點做坐標轉換，再重新映射紋理，最后保存模型到文件。但該方法非常復雜，實現難度大。最簡單直觀的思路是直接將OSGB文件中的投影坐標替換為轉換并偏移后的坐標，并修改metadata.xml和S3C文件中的坐標系及坐標偏移量信息。

2.1 待轉換坐標獲取

在Geometry類對象的屬性中，一般只有頂點數組含有需要轉換的坐標數據。PagedLOD類對象根據視點到物體中心的距離決定顯示哪種精細度的模型，而確定物體中心的方法由屬性“中心模式”給出，既可以是包圍球中心，也可以是用戶自定義中心，還可以是二者的結合。對于后兩種模式，需要對用戶自定義中心的坐標進行轉換，否則可能會造成模型顯示時被誤裁剪。文中定義三維模型需要轉換的坐標包括Geometry類對象的頂點坐標和PagedLOD類對象的用戶自定義中心坐標。

(1)

2.2 投影坐標轉換

將投影坐標轉換分為相同橢球體的投影坐標轉換和不同橢球體間的投影坐標轉換。前者的轉換框架可見圖4。后者采用七參數法[9]進行轉換，其中未假設3個旋轉角很小，轉換框架可見圖5。轉換前后高程坐標不變。

圖4 相同橢球體的投影坐標轉換框架

圖5 不同橢球體間的投影坐標轉換框架

2.3 轉換模型保存

(2)

在坐標轉換和偏移完成后，需要將轉換后模型保存到OSGB文件。理論上，可以調用osgDB::Registry類實例的writeNode函數將模型保存為OSGB文件。但這樣做存在格式兼容性風險。如果ContextCapture不能支持所有版本的OSGB格式，則該方法生成的OSGB文件可能不能被ContextCapture準確讀取。

一種簡單的方法是直接對原始OSGB文件進行坐標數據替換。在OSGB文件的數據流中，一個Geometry類對象的頂點坐標連續放置，可被視為一個浮點型數組，而PagedLOD類對象的用戶自定義中心坐標可被視為一個具有3個元素的浮點型數組。設轉換前的數組為Asrc，存儲了N個三維坐標，元素數目為3N，對應的轉換和偏移后的數組為Adst，則以字節為單位的Asrc和Adst的大小nLength均為12N(單精度浮點型)或24N(雙精度浮點型)。通過在文件數據流中比較長度為nLength的連續內存的值是否和Asrc相同來確定數據替換位置，在找到替換位置后，用Adst的值覆蓋原始數據。

3 實現方法

利用C++和開源軟件開發OSGB模型的投影坐標轉換軟件，其中Qt[13]用作編程框架，OpenSceneGraph 3.7[4]用于三維模型讀取，PROJ 6.0[14]用于地圖投影，GDAL 3.0[15]用于坐標系設置，Eigen3[16]用于矩陣運算。

3.1 坐標轉換流程

圖 6為設計的坐標轉換流程。

圖6 模型投影坐標轉換流程

3.2 轉換前模型解析

3.2.1 源坐標系解析

如果模型源坐標系以EPSG編碼的形式提供，可調用OGRSpatialReference類的importFromEPSG函數導入編碼對應的坐標系；如果模型源坐標系以WKT文本的形式提供，可調用OGRSpatialReference類的importFromWkt函數導入WKT對應的坐標系。

3.2.2 OSGB文件解析

要讀取OSGB文件中的模型，可先調用osgDB::readNodeFile函數獲得Node類型的模型對象。之后，要獲取待轉換坐標，需要對OSG定義的類NodeVisitor進行派生并重寫相應的apply虛函數。示例代碼如下：

class CTVisitor: public NodeVisitor {

public:

virtual void apply(Geode &node);

virtual void apply(Group &node);

}

對Geometry類對象的解析可放在void CTVisitor::apply(Geode &node)中，示例如下：

for(int n = 0; n < node.getNumDrawables(); n++){

Drawable *pDb = node.getDrawable(n);

if(pDb->className()== string("Geometry")){

ref_ptr pGm =(Geometry*)pDb;

//pVt為頂點坐標數組指針

Vec3Array *pVt =dynamic_cast

(pGm->getVertexArray().get());

}

}

traverse(node);

如果模型采用雙精度浮點數存儲頂點坐標，則需要將上述代碼中的Vec3Array替換為Vec3dArray。對PagedLOD類對象的解析可放在void CTVisitor::apply(Group &node)中，示例如下：

if(node.className()== string("PagedLOD")){

PagedLOD lod =(PagedLOD)node;

if((lod.getCenterMode()== LOD::USER_ DEFINED_CENTER)||(lod.getCenterMode()== LOD::UNION_OF_BOUNDING_SPHERE_AND _USER_DEFINED)){

//cent為用戶自定義中心坐標

Vec3 cent = lod.getCenter();

}

}

traverse(node);

3.3 投影坐標轉換

3.3.1 大地坐標與投影坐標間的轉換

大地坐標與投影坐標間的轉換采用pj_transform函數，該函數定義如下：

int pj_transform(projPJ src, projPJ dst, long pt_count, int pt_offset, double *x, double *y, double *z)

其中，src為源坐標系對象，dst為目標坐標系對象，pt_count為待轉換點數，pt_offset為首個待轉換點在數組中的序號，x，y，z分別為存儲三維坐標的數組的地址。在PROJ中，大地和投影坐標系對象需用PROJ字符串定義。大地坐標系的PROJ字符串為：

QString("+proj=longlat +ellps=%1 +no_defs +type= crs").arg(strEllps)

高斯-克呂格投影坐標系的PROJ字符串為：

QString("+proj=tmerc +lat_0=0 +lon_0=%1 +k=1 +x_0=%2 +y_0=%3 +ellps=%4 +units=m +no_defs +type=crs").arg(dCentLon).arg(dFalseEast).arg(dFalseNorth).arg(strEllps)

其中，strEllps為大地坐標系使用的橢球體的字符串，1954年北京坐標系是“krass”，1980年西安坐標系是“IAU76”，2000國家大地坐標系是“GRS80”；dCentLon為中央子午線以角度制表示的經度，dFalseEast為以m為單位的東偏，dFalseNorth為以m為單位的北偏。之后，調用pj_init_plus函數，輸入參數為坐標系的PROJ字符串，即可獲得projPJ類型的坐標系對象。

3.3.2 大地坐標與空間直角坐標間的轉換

從大地坐標轉換到空間直角坐標采用pj_geodetic_to_geocentric函數，該函數定義如下：

int pj_geodetic_to_geocentric(double a, double es, long pt_count, int pt_offset, double *x, double *y, double *z)

從空間直角坐標轉換到大地坐標采用pj_geocentric_to_geodetic函數，該函數定義如下：

int pj_geocentric_to_geodetic(double a, double es, long pt_count, int pt_offset, double *x, double *y, double *z)

其中，a為橢球體半長軸，es為第一偏心率平方，pt_count為待轉換點數，pt_offset為首個待轉換點在數組中的序號，x，y，z分別為存儲三維坐標的數組的地址。

3.4 轉換后模型生成

在生成轉換后模型時，需要確保轉換前后目錄結構不變，OSGB文件名和相對路徑不變。

3.4.1 metadata.xml生成

在轉換后模型的metadata.xml中，對于EPSG已定義的投影坐標系，可將SRS標簽值設置為坐標系的EPSG編碼；對于EPSG未定義的投影坐標系，可先調用OGRSpatialReference類的SetWellKnownGeogCS函數設置大地坐標系，再調用SetTM函數設置中央子午線經度、中央緯線緯度、中央子午線縮放因子、東偏、北偏等參數，調用SetProjCS函數設置坐標系名稱，調用exportToWkt函數獲得坐標系的WKT，最后將SRS標簽值設置為WKT。SRSOrigin標簽值為Δxdst，Δydst,Δzdst。

3.4.2 OSGB文件生成

設Asrc的數組指針為pSrc，Adst的數組指針為pDst，以字節為單位的OSGB文件大小為nOsgbSize。先將轉換前的OSGB文件全部讀入unsigned char類型、大小為nOsgbSize的數組中。設數組指針為pszOsgb，則坐標數據替換的示例代碼如下：

for(int n = 0; n < nOsgbSize-nLength; n++){

int b = memcmp(pszOsgb + n, pSrc, nLength);

if(b == 0){

memcpy(pszOsgb + n, pDst, nLength);

break;

}

}

最后將pszOsgb內的數據保存到轉換后目錄下的OSGB文件。

3.4.3 S3C文件坐標系設置

雖然S3C文件也存儲坐標系信息，但由于新版S3C文件的格式已做加密，難以用編程的方法修改其坐標系信息。因此，先將轉換前的S3C文件復制到轉換后目錄下，再利用CC_S3CComposer.exe將S3C文件中的坐標系信息替換為metadata.xml中的坐標系信息。

4 軟件實例

圖7是軟件界面，支持中國常用大地坐標系基礎上的高斯-克呂格投影；支持自定義和自動計算輸出坐標偏移量；支持手動和文件輸入橢球體轉換七參數。

圖7 軟件界面

在Intel Core i5-2300的CPU、24GB內存的臺式機上，軟件將具有6個瓦片、1179個OSGB文件的三維模型從2000國家大地坐標系基礎上120°中央子午線的投影坐標系轉換到1954年北京大地坐標系基礎上120.75°中央子午線的投影坐標系，耗時20.214 s，平均單個OSGB文件耗時17 ms。在利用CC_S3CComposer.exe更新S3C文件的坐標系信息后，可以用ContextCapture Viewer[1]打開轉換后的模型，幾何無變形，紋理無缺失，顯示無閃爍。為了對比，由人工在目標坐標系下重新生產三維模型。通過在軟件轉換后模型和人工重新生產的模型上選取同名特征點并查詢其坐標，發現坐標一致。

5 結束語

文中提出的投影坐標轉換方法簡單高效，對OSGB文件的格式要求不高，不受是否采用LOD技術、有無紋理貼圖、紋理壓縮比例、切塊方式等因素的影響。由于OSG基本支持所有主流的三維模型格式，所以文中方法容易擴展到其它二進制的三維模型格式。此外，坐標數據替換法也適用于其它需要對模型坐標進行修改的情況。

文中方法未修改Geometry類對象的法向量，但理論上，修改頂點坐標之后需要重新計算法向量。雖然在第4節例子中不更新法向量在視覺上沒有明顯影響，但更新模型法向量仍是一個改進方向。