1 引言

　　隨著天文觀測設(shè)備的不斷發(fā)展，人類對宇宙的探索能力也不斷提高。特別是最近幾十年來，計(jì)算機(jī)和電荷耦合器件 (charge-coupled device, CCD) 在天文學(xué)領(lǐng)域的廣泛使用使得人類從過去的肉眼觀測、膠片觀測進(jìn)入到數(shù)字化信息時(shí)代。多個(gè)大規(guī)模數(shù)字巡天項(xiàng)目預(yù)期能夠獲得的數(shù)據(jù)量遠(yuǎn)遠(yuǎn)超過人類手工處理的極限。例如建設(shè)在智利的 8.4m 口徑望遠(yuǎn)鏡薇拉・魯賓天文臺(tái) (Vera C. Rubin Observatory Legacy Survey of Space and Time, LSST) 計(jì)劃從 2022 年開始開展為期 10a 的南半球 18000 平方度巡天。在 10a 中，每一個(gè)天區(qū)將會(huì)被訪問 1000 次以上，獲得 400 億天體的觀測數(shù)據(jù)以及一個(gè)空前規(guī)模的動(dòng)態(tài)時(shí)域巡天數(shù)據(jù)庫。泛星計(jì)劃 (Panoramic Survey Telescope and Rapid Response System, Pan-STARRS) 使用 4 個(gè)口徑為 1.8m 的望遠(yuǎn)鏡組成陣列。

　　一方面，這些大規(guī)模巡天觀測將產(chǎn)生數(shù)以太字節(jié) (terabyte, TB) 乃至拍字節(jié) (petabyte, PB) 量級(jí)的數(shù)據(jù)。例如，Pan-STARRS 基于第一臺(tái)望遠(yuǎn)鏡 Pan-STARRS1 (PS1)，該項(xiàng)目已發(fā)布 DR1、DR2 兩批數(shù)據(jù)。在 2019 年 1 月 28 日發(fā)布的 PS1 DR2 星表中，ObjectThin 表包含 105 億行數(shù)據(jù) (約 5.4TB)，StackObjectThin 表包含 34.7 億行數(shù)據(jù) (約 1.2TB)。這些數(shù)據(jù)需要及時(shí)處理與分析。另一方面，其中的暫現(xiàn)源搜尋類巡天每次觀測得到的數(shù)據(jù)又需要及時(shí)與之前大量的觀測數(shù)據(jù)交叉證認(rèn)，以判斷是否有變化的天體。

　　前人在星表檢索方面的很多工作都基于數(shù)據(jù)庫而實(shí)現(xiàn)。對于天體坐標(biāo)這類空間二維數(shù)據(jù)而言，若直接將天體信息存入數(shù)據(jù)庫，在使用結(jié)構(gòu)化查詢語言 (structured query language, SQL) 進(jìn)行錐形檢索時(shí)，需要通過復(fù)雜的三角函數(shù)計(jì)算各個(gè)坐標(biāo)之間的角距離。少部分?jǐn)?shù)據(jù)庫如 PostgreSQL 有 PostGIS 等擴(kuò)展組件支持對空間數(shù)據(jù)的檢索，但是這些組件為地理信息系統(tǒng) (geographic information system, GIS) 所設(shè)計(jì)，更多地用于處理地理上的三維點(diǎn)線面等復(fù)雜區(qū)域，而天文星表中的天體多為二維點(diǎn)源，且涉及天文坐標(biāo)與地理坐標(biāo)的轉(zhuǎn)換等操作，在天文上少見應(yīng)用。

　　更普遍的做法是使用偽二維球面索引算法將天球劃分成多個(gè)子區(qū)域，然后將每個(gè)區(qū)域賦予一個(gè) id 編號(hào)，星表中的天體以所在區(qū)域的 id 為主鍵存入數(shù)據(jù)庫中。在檢索時(shí)根據(jù)坐標(biāo)計(jì)算區(qū)域 id 讀取該區(qū)域的天體，再計(jì)算距離。目前球面索引算法有很多，例如條帶算法 (Zones Algorithm)、分層三角網(wǎng)格 (Hierarchical Triangular Mesh, HTM)、四叉樹立方體 (Quad Tree Cube, Q3C)、多級(jí)等面積同緯度劃分法 (Hierarchical Equal Area isoLatitude Pixelisation, HEALPix) 等。Berriman 等人使用 2MASS 全天源星表 (two micron all sky survey all-sky point source catalog )(約 4.7 億天體) 和未合并的 HSC 星表 (non-merged Hubble source catalog)(約 3.83 億天體) 測試了 HTM、HEALPix 的不同層級(jí)分別在 Solaris、Windows、Red Hat Linux 平臺(tái)下的 PostgreSQL 數(shù)據(jù)庫和 SQL Server 數(shù)據(jù)庫的檢索速度。Szalay 等人開發(fā)的 SkyServer 采用微軟的 SQL Server 數(shù)據(jù)庫和 HTM 算法，實(shí)現(xiàn)了對斯隆數(shù)字巡天 (Sloan Digital Sky Survey, SDSS) 數(shù)據(jù)的在線檢索和訪問。這也是第一次使用商業(yè)數(shù)據(jù)庫實(shí)現(xiàn)了大規(guī)模天文數(shù)據(jù)管理與訪問。

　　星表檢索的一個(gè)重要應(yīng)用就是實(shí)現(xiàn)對多個(gè)星表的交叉證認(rèn)與融合，高丹通過使用 HTM 和數(shù)據(jù)庫實(shí)現(xiàn)了對數(shù)十萬量級(jí)星表的檢索與交叉證認(rèn)，但沒能解決劃分邊界附近天體的漏源問題。Du 等人通過同時(shí)使用 HEALPix 和 HTM 以及數(shù)據(jù)庫實(shí)現(xiàn)了對星表的交叉證認(rèn)，緩解了漏源問題。Boch 等人通過使用 KD-Tree 和 HEALPix 實(shí)現(xiàn)了在線交叉證認(rèn)平臺(tái) CDS-Xmatch，該平臺(tái)允許用戶在網(wǎng)頁上傳自定義星表與平臺(tái)自帶或之前上傳的星表進(jìn)行交叉證認(rèn)。

　　前人的工作已經(jīng)證實(shí)使用數(shù)據(jù)庫可以實(shí)現(xiàn)較高效的星表檢索。在使用偽二維球面索引算法對星表進(jìn)行基于坐標(biāo)的檢索時(shí)，數(shù)據(jù)庫的作用是存儲(chǔ)天體信息和對應(yīng) id 編號(hào)，在檢索時(shí)需要找到所有與待檢索 id 相同的天體，交由上層軟件計(jì)算角距離并加以過濾后，才能得到最終符合要求的結(jié)果。根據(jù)層級(jí)選擇的不同，可能會(huì)有成千上萬個(gè)天體都具有相同的區(qū)域 id，在這種特定場合下，若不使用數(shù)據(jù)庫，而是直接將相同 id 的所有天體保存在同一文件中，檢索時(shí)只需讀取對應(yīng)的文件即可得到該 id 編號(hào)下所有的天體，也可以實(shí)現(xiàn)對星表的索引與檢索功能。

　　在劃分時(shí)區(qū)域編號(hào)對于數(shù)據(jù)庫僅是一個(gè)用于建立索引的數(shù)字，但是 HEALPix 等類似算法的不同層級(jí)、同一層級(jí)的不同區(qū)域之間的 id 編號(hào)往往有一定的關(guān)聯(lián)，在使用文件作為星表的存儲(chǔ)介質(zhì)時(shí)，可以利用這些特點(diǎn)做一些針對星表的優(yōu)化，提升索引與檢索的效率。此外，與數(shù)據(jù)庫這類商業(yè)軟件不同，文件具有無需授權(quán)，通用性更好，也無需提前安裝配置，易于維護(hù)等優(yōu)點(diǎn)。在諸如偏遠(yuǎn)臺(tái)站等較難維護(hù)，或者在單片機(jī)、小型探測器等這類性能較低，難以運(yùn)行數(shù)據(jù)庫的場合，文件都具有獨(dú)特的優(yōu)勢。因此，本文嘗試使用文件作為星表檢索時(shí)的存儲(chǔ)介質(zhì)，并加以優(yōu)化，以提高在星表交叉證認(rèn)或融合時(shí)的檢索效率。

　　本文選擇了 HEALPix 作為劃分的基礎(chǔ)，該方法由 Górski 等人提出，能夠根據(jù)預(yù)先設(shè)定的層級(jí)將天球切分成多個(gè)面積相同的子區(qū)域，稱之為 “像素”，每個(gè)像素都有對應(yīng)的編號(hào)。層級(jí)越大，像素越多，單位像素面積越小，對天球劃分的分辨率越高。HEALPix 最開始應(yīng)用于宇宙微波背景輻射的研究與分析，如今已經(jīng)在天文學(xué)研究中廣泛使用。Pineau 等人使用 HEALPix 作為劃分方案對星表進(jìn)行了交叉證認(rèn)測試。原始 HEALPix 只能以固定分辨率劃分球面，之后有各種多分辨率的 HEALPix 解決方案被提出，例如多層級(jí)覆蓋天區(qū) (multi-order coverage map, MOC)、多分辨率 HEALPix (multi-resolution healpix, MRH)、mhealpy 等。許允飛通過使用 MOC 建立的多層級(jí)覆蓋天區(qū)四叉樹 (multi-order coverage tree, MOC-Tree)，從而對不規(guī)則區(qū)域星表進(jìn)行索引。Martinez-Castellanos 等人發(fā)布的 mhealpy 是一套多分辨率 HEALPix 及可視化的 Python 包，支持兩種多分辨率模式，其中一種允許用戶通過自定義函數(shù)選擇當(dāng)前區(qū)域的劃分層級(jí)，但是受到可視化的限制，某個(gè)區(qū)域只能被一個(gè)層級(jí)所覆蓋。而星表的切分不需要考慮可視化相關(guān)的限制，本文嘗試允許多個(gè)層級(jí)的像素覆蓋同一區(qū)域，實(shí)現(xiàn)對星表的均勻切分。

　　天文星表在計(jì)算機(jī)中存儲(chǔ)和交換的格式以逗號(hào)分隔值 (comma-separated values, CSV) 和靈活圖像傳輸系統(tǒng) (flexible image transport system, FITS) 為主，此外常用的還有 VOTable 和一些機(jī)構(gòu)制定的文本或二進(jìn)制格式。CSV 格式的星表以純文本格式存儲(chǔ)，使用逗號(hào)或其他符號(hào)進(jìn)行分隔，使用文本編輯器即可讀取和修改內(nèi)容。程序使用 CSV 存取星表的天體坐標(biāo)時(shí)需要進(jìn)行字符串和浮點(diǎn)數(shù)的轉(zhuǎn)換，對存取大規(guī)模星表的速度有一定的影響。

　　FITS 是由國際天文學(xué)聯(lián)合會(huì)定義的信息交換格式，最初是為了傳輸圖像而設(shè)計(jì)，之后擴(kuò)展到星表等更復(fù)雜的數(shù)據(jù)格式;FITS 格式的星表以二進(jìn)制格式按列存儲(chǔ)，讀寫需要依靠第三方軟件或庫。

　　在檢索星表時(shí)僅需要讀取坐標(biāo)等少量數(shù)據(jù)，也不需要考慮可讀性。另外這些文件僅供檢索程序運(yùn)行時(shí)使用，不需要對外分享。可以考慮通過定制專屬的序列化協(xié)議，對這些數(shù)據(jù)序列化，以減少中間格式轉(zhuǎn)換的次數(shù)，達(dá)到更好的存取性能。

　　巡天觀測掃描的天空區(qū)域通常是連續(xù)的，前一次檢索的部分?jǐn)?shù)據(jù)可為下一次檢索所用。通過選擇適當(dāng)?shù)木彺娌呗裕梢詼p少文件讀取次數(shù)，提高數(shù)據(jù)的復(fù)用率，并實(shí)現(xiàn)更快的連續(xù)檢索。除此之外，在執(zhí)行多個(gè)大規(guī)模星表融合等需要對大量天體交叉證認(rèn)的任務(wù)時(shí)，通過選取合適的遍歷順序，可以提高緩存的命中率，提高交叉證認(rèn)的速度。

　　針對在基于文件的情況下如何對星表高效檢索的問題，本文提出了一套解決方案。本文第 2 章研究如何在有限內(nèi)存下將大規(guī)模星表切分成合適且均勻的星表子集;第 3 章研究如何對切分后的星表子集做高效檢索;第 4 章探討如何使用開源序列化庫對星表高效存取;第 5 章將討論實(shí)現(xiàn)多個(gè)大規(guī)模星表的交叉證認(rèn)與融合的方法;第 6 章將對前幾章提出的方法做實(shí)驗(yàn)測試;第 7 章進(jìn)行總結(jié)。

　　2 星表索引及星表文件的切分策略

　　從一個(gè)包含有 X 個(gè)天體的星表中搜尋 Y 個(gè)特定天體，若不對星表做任何處理，直接從頭到尾按順序查找，最壞情況下需 X×Y 次檢索。一個(gè)改進(jìn)方法是用諸如 KD-Tree 等多維數(shù)據(jù)結(jié)構(gòu)建立索引加快檢索時(shí)的速度。受制于當(dāng)前計(jì)算機(jī)體系結(jié)構(gòu)，這些數(shù)據(jù)結(jié)構(gòu)建立和檢索只能在內(nèi)存中進(jìn)行。但對于動(dòng)輒數(shù) TB 的星表，一方面計(jì)算機(jī)內(nèi)存難以將星表全部讀入到內(nèi)存，另一方面從硬盤中將星表完整讀取一遍也會(huì)消耗大量時(shí)間。為了加快查找速度，需要針對性的優(yōu)化。

　　對于星表來說，有以下兩種基本的基于坐標(biāo)的查詢操作：

　　最近鄰天體檢索，即給定坐標(biāo)，檢索星表中與之最為接近的天體;

　　錐形檢索，即給定坐標(biāo)和距離閾值，檢索星表中所有距離該坐標(biāo)小于該閾值的天體子集。

　　這兩種操作只需要檢索給定坐標(biāo)附近區(qū)域內(nèi)天體，并不需要將整個(gè)星表都讀取并檢索一遍。

　　因此，在檢索前可以將星表按照天體的坐標(biāo)執(zhí)行切分，把鄰近的天體切分到同一個(gè)文件中，并給切分到的每個(gè)文件賦予一個(gè)特定的編號(hào)。在檢索時(shí)只需要根據(jù)待檢索坐標(biāo)等信息，計(jì)算所覆蓋到區(qū)域的編號(hào)，并根據(jù)該編號(hào)讀取并檢索相應(yīng)的文件，這樣就可以避免讀取整個(gè)星表，從而大大提高查詢效率。

　　切分出的星表子集應(yīng)做到大小合適且盡量均勻。若單個(gè)文件內(nèi)天體數(shù)量太多，則檢索時(shí)讀取時(shí)間和建立索引的時(shí)間較長，影響檢索的效率;若數(shù)量太少，則每次檢索需要加載大量的星表文件，而且硬盤對小文件的隨機(jī)讀寫速度遠(yuǎn)遠(yuǎn)小于順序讀寫，檢索時(shí)需要跨越多個(gè)星表文件，極限情況下甚至?xí)嘶身樞驒z索，造成效率的下降。

　　因而，需要對切分后單個(gè)星表文件天體數(shù)量上下限做出限制，以最大化星表檢索效率。

　　2.1 多分辨率 HEALPix 星表切分算法

　　在星表中，天體的位置坐標(biāo)作為必不可少的信息之一，適合作為切分的首選依據(jù)。但直接按照赤經(jīng)和赤緯的數(shù)值線性劃分得到的區(qū)域面積不等，難以實(shí)現(xiàn)均勻的劃分，本文選擇能將天球等面積劃分的 HEALPix 作為切分的基礎(chǔ)。

　　雖然 HEALPix 對天球的劃分是等面積的，但是天體在天球中的分布往往是不均勻的。例如，對于一個(gè)以河內(nèi)恒星為主的星表來說，銀心方向的天體密度往往大于反銀心方向或者垂直于銀道面方向。原始的 HEALPix 只能以一個(gè)固定的層級(jí)劃分天球，導(dǎo)致不同區(qū)域天體數(shù)量相差較大，以此切分得到的星表文件也會(huì)有較大的體積差異，影響對星表處理的效率。為了保證切分出來的星表盡可能均勻，本文提出了多個(gè) HEALPix 層級(jí)共存的星表切分算法，能夠根據(jù)區(qū)域天體的密度動(dòng)態(tài)選取一個(gè)合適的切分層次。

　　HEALPix 首先將天球劃分成 12 個(gè)等面積的大區(qū)域，并在此基礎(chǔ)上根據(jù)一個(gè)設(shè)定的層級(jí) (設(shè)為 k 將整個(gè)天球劃分成 12×4ᵏ個(gè)等面積的小區(qū)域。k 值越大，劃分精度越細(xì)，分辨率越高。

　　HEALPix 支持 RING 和 NESTED 兩種編號(hào)模式。RING 模式下，從北極到南極，相同緯度的區(qū)域形成一個(gè)環(huán)，在環(huán)內(nèi)按照順序依次編號(hào)。NESTED 模式下，將第 i 層編號(hào)為 nₚᵢₓ,ᵢ的區(qū)域繼續(xù)劃分成 4 份，即可得到 4 個(gè)第 i+1 層的區(qū)域，編號(hào)范圍是 [4×nₚᵢₓ,ᵢ,4×nₚᵢₓ,ᵢ+3]。這樣的特性非常適合應(yīng)用四叉樹。因此，本文選擇在 NESTED 模式下使用四叉樹作為聯(lián)系起不同層次區(qū)域之間的索引方案。

　　四叉樹是一種基于空間遞歸分解思想的分層數(shù)據(jù)結(jié)構(gòu)。樹的根節(jié)點(diǎn)可以代表總的區(qū)域，將該區(qū)域劃分成 4 個(gè)子區(qū)域后，每個(gè)子區(qū)域都可以用四叉樹的一個(gè)子節(jié)點(diǎn)來表示。每個(gè)子區(qū)域還可以再細(xì)分為四個(gè)更小的區(qū)域，這個(gè)過程可以一直遞歸進(jìn)行，直到達(dá)到所需的分辨率。四叉樹與 HEALPix 層級(jí)之間的映射關(guān)系如圖 3 所示。

　　由于 HEALPix 不同層級(jí)區(qū)域都是由最初的 12 個(gè)區(qū)域劃分而來，所以可建立 12 棵四叉樹，每棵樹的根節(jié)點(diǎn)代表 k=0 時(shí)的 12 個(gè)區(qū)域其中之一。樹第 i 層的節(jié)點(diǎn)，與 k=i 時(shí)的劃分區(qū)域一一對應(yīng)。樹的每個(gè)非葉節(jié)點(diǎn)都有 4 個(gè)子節(jié)點(diǎn)，指向下一層的 4 個(gè)子天區(qū)。每個(gè)子節(jié)點(diǎn)的 nₚᵢₓ編號(hào)與父節(jié)點(diǎn)有相同的二進(jìn)制前綴。這樣就實(shí)現(xiàn)了四叉樹與 HEALPix 多個(gè)層級(jí)劃分區(qū)域的一一對應(yīng)，可以使用對樹的相關(guān)操作來實(shí)現(xiàn)對天區(qū)的有序訪問。

　　由于 HEALPix 在不同層級(jí)下對區(qū)域均是從 0 開始編號(hào)，因此不同層級(jí)之間的區(qū)域編號(hào)會(huì)產(chǎn)生沖突，這種二維編號(hào)方式會(huì)影響檢索的效率。本文采用 MOC 中的 NUNIQ 編碼方式，賦予不同層級(jí)的區(qū)域一個(gè)唯一的編號(hào)。

　　由于大規(guī)模星表大小 N 往往會(huì)遠(yuǎn)遠(yuǎn)大于可用內(nèi)存容量 M，而程序運(yùn)行時(shí)所能讀取的星表又不能大于 M。此時(shí)就需要將星表事先分割成⌈N/M⌉塊局部星表，每塊大小都小于 M。這樣，程序每次都可以將一塊局部星表完整讀入到內(nèi)存執(zhí)行切分，之后清空內(nèi)存并讀取下一塊局部星表。而對多個(gè)局部星表的切分后直接的合并，并不等同于對原始星表的切分結(jié)果，需要再次利用四叉樹歸并。所以本文的星表切分算法分成兩步：第一步是對星表的局部切分，第二步是將多個(gè)局部切分的星表歸并。

　　2.2 星表的局部切分

　　首先需要指定一個(gè)參數(shù) kₘₐₓ，代表所能使用的 HEALPix 的最大層級(jí)，并作為四叉樹的最大層級(jí)。但 kₘₐₓ越大，程序運(yùn)行時(shí)消耗的內(nèi)存也越多。

　　第一步是四叉樹的建立與星表的讀取。如上文所述，建立 12 棵四叉樹，每棵樹有 kₘₐₓ層，第 i (i∈[0,kₘₐₓ-1]) 層有 4ⁱ個(gè)節(jié)點(diǎn)。12 棵樹共計(jì) 4ᵏ⁽ᵐᵃˣ⁾⁺¹-4 個(gè)節(jié)點(diǎn)，每個(gè)節(jié)點(diǎn)都有一個(gè) UNIQ 編號(hào)并與天球上一個(gè)區(qū)域相對應(yīng)。

　　四叉樹的每個(gè)節(jié)點(diǎn)內(nèi)有一個(gè)變量 C 代表該節(jié)點(diǎn)對應(yīng)區(qū)域內(nèi)的天體數(shù)量，初始時(shí)為 0。每讀取星表中一個(gè)天體信息，就根據(jù)它的坐標(biāo)值，分別計(jì)算出該坐標(biāo)在每個(gè)層級(jí)下的 HEALPix 的區(qū)域編號(hào)。并根據(jù)此編號(hào)計(jì)算對應(yīng)的 UNIQ 編號(hào)，同時(shí)將對應(yīng)的節(jié)點(diǎn) c 值加 1。即實(shí)現(xiàn)了天體在四叉樹中的索引。

　　當(dāng)程序?qū)⒁粔K局部星表讀入到內(nèi)存并用四叉樹索引后，將執(zhí)行切分操作。切分算法需要事先指定一個(gè)切分閾值 T，在不同切分模式下作為單個(gè)星表文件的天體數(shù)量上限或者下限，當(dāng)遍歷到某一節(jié)點(diǎn)時(shí)，根據(jù)該節(jié)點(diǎn)所覆蓋區(qū)域內(nèi)的天體總數(shù)與閾值 T 相比較，以確定下一步的操作。在遞歸切分時(shí)，會(huì)遇到以下 3 種情況：

　　父節(jié)點(diǎn)的 c 值和 4 個(gè)子節(jié)點(diǎn)中的 C 值均大于 T，則繼續(xù)遞歸遍歷 4 個(gè)子節(jié)點(diǎn);

　　父節(jié)點(diǎn)的 c 值大于 T，且存在子節(jié)點(diǎn)的 c 值小于 T 的節(jié)點(diǎn)，處理方式見下文;

　　已到達(dá) kₘₐₓ層，則將當(dāng)前節(jié)點(diǎn)內(nèi)的天體切分，并停止向下遍歷。

　　在星表分布不均勻的情況下，例如北半球地面望遠(yuǎn)鏡受地理緯度的限制難以觀測到南天區(qū)的天體，會(huì)出現(xiàn)南天區(qū)中根節(jié)點(diǎn) C 值也小于 T 的情況，此時(shí)為特例情況，將根節(jié)點(diǎn)內(nèi)所有天體切分即可。

　　對于情況 (2)，又可以分為以下 3 種處理方式，如圖 5 舉例所示：

　　“激進(jìn)” 方式。對 c 值大于 T 的子節(jié)點(diǎn)繼續(xù)遍歷其子節(jié)點(diǎn);對小于 T 的子節(jié)點(diǎn)，則將子節(jié)點(diǎn)對應(yīng)區(qū)域的天體切分至對應(yīng)文件中，并停止該子節(jié)點(diǎn)的向下遞歸，即在圖 5 例中將 4Xp 和 4×p+3 的節(jié)點(diǎn)各自切分到兩個(gè)文件中。該方式是 mhealpy 中所提供的 2 個(gè)切分方案之一，在本文中按照其特點(diǎn)將其命名為 “激進(jìn)” 方式;

　　“保守” 方式。只要 4 個(gè)子節(jié)點(diǎn)中存在 1 個(gè) c 值小于 T 的節(jié)點(diǎn)，就將父節(jié)點(diǎn)對應(yīng)區(qū)域的天體切分到對應(yīng)文件中，并停止遍歷其子節(jié)點(diǎn)，即將圖 5 中編號(hào)為 p 的父節(jié)點(diǎn)切分;

　　“中立” 方式。分別訪問該節(jié)點(diǎn)的 4 個(gè)子節(jié)點(diǎn)，若該節(jié)點(diǎn) c 值大于 T，則繼續(xù)遍歷其子節(jié)點(diǎn)，并得到一個(gè)返回值 r，代表了對應(yīng)區(qū)域有 r 個(gè)天體被切分。若該節(jié)點(diǎn) C 值小于 T，則不再遍歷其子節(jié)點(diǎn)，而是將 C 的值作為 r 返回至父節(jié)點(diǎn)。父節(jié)點(diǎn)收到各個(gè)子節(jié)點(diǎn)的返回值后計(jì)算剩余未被切分的天體數(shù)目，若總數(shù)大于 T，則以父節(jié)點(diǎn)的名義將這些天體切分到對應(yīng)文件中，并返回總數(shù);若小于 T，則將已切分的數(shù)量返回至更上一層節(jié)點(diǎn)。依次類推，直到該節(jié)點(diǎn)總數(shù)大于 T，或到達(dá)根節(jié)點(diǎn)時(shí)將天體切分到文件中。圖中假設(shè) 4×p+1 和 4×p+2 區(qū)域已全被切分。

　　經(jīng)過對比不難看出，本文所提出的 “中立” 切分方式在理論上可以切分出最均衡的星表文件，適合作為實(shí)際應(yīng)用時(shí)的首選切分方式。在切分完成之后，清空內(nèi)存中的星表，若星表未讀取完畢，則繼續(xù)讀取星表，并重復(fù)上述切分步驟，直至所有星表讀取完畢，完成星表的局部切分。

　　2.3 局部切分星表的歸并

　　一個(gè)無序的星表經(jīng)過⌈N/M⌉次局部切分后，得到了⌈N/M⌉塊局部切分的星表文件集合，每塊集合內(nèi)有不定數(shù)量的星表文件，雖然塊內(nèi)的天體都按照所在區(qū)域劃分到各自的文件中，但分塊與分塊之間卻是互相獨(dú)立的。因此需要?dú)w并各個(gè)分塊的星表，得到對原始星表完整的切分結(jié)果。

　　同樣是建立 12 棵對應(yīng)的四叉樹。從每棵樹的根節(jié)點(diǎn)開始，按照深度優(yōu)先的順序遍歷樹的各個(gè)子節(jié)點(diǎn)。若某個(gè)局部切分的星表中存在該節(jié)點(diǎn)對應(yīng)的星表文件，則讀取到內(nèi)存，并更新對應(yīng)四叉樹節(jié)點(diǎn)中的 c 值。之后依次遞歸遍歷該節(jié)點(diǎn)的 4 個(gè)子節(jié)點(diǎn)，并接受各個(gè)子節(jié)點(diǎn)的返回值(僅 “中立” 模式有返回值)。

　　當(dāng)遍歷完 4 個(gè)子節(jié)點(diǎn)并返回到當(dāng)前節(jié)點(diǎn)時(shí)，則意味著以該節(jié)點(diǎn)為根的所有子節(jié)點(diǎn)均已被訪問過，星表中該節(jié)點(diǎn)所在區(qū)域的所有天體均已被訪問過，或留在內(nèi)存中，或已經(jīng)在遍歷某個(gè)子節(jié)點(diǎn)時(shí)被切分到最終星表文件中。此時(shí)，與局部切分過程類似，根據(jù)切分模式、閾值 T 和當(dāng)前節(jié)點(diǎn)所剩余的天體數(shù)量，判斷是否執(zhí)行切分。“中立” 模式還需向父節(jié)點(diǎn)返回對應(yīng)的值。此時(shí)切分得到的星表文件即為歸并之后的星表文件。此過程不斷進(jìn)行，最后實(shí)現(xiàn)對局部星表的歸并。

　　可以看到，天體密集的區(qū)域切分的精細(xì)程度比較高，天體稀疏的區(qū)域切分精細(xì)度較低，保持了每個(gè)切分區(qū)域天體數(shù)量的均衡。

　　3 切分后星表的檢索

　　在將星表按照上述切分算法切分之后，每個(gè)切分后的星表文件都有一個(gè) UNIQ 編號(hào)，代表了其所覆蓋天區(qū)。在檢索時(shí)，根據(jù)不同檢索方式，將待檢索區(qū)域轉(zhuǎn)換為對應(yīng)的 UNIQ 編號(hào)列表并讀取對應(yīng)的文件。

　　3.1 區(qū)域檢索

　　由于 “中立” 模式下多重覆蓋的特性，對于最高層級(jí)的某個(gè)區(qū)域來說，該區(qū)域內(nèi)的天體只會(huì)存在于該區(qū)域所在的節(jié)點(diǎn)或距離該節(jié)點(diǎn)最近的有星表子集存在的上層節(jié)點(diǎn)中。而一個(gè)非最高層級(jí)的區(qū)域內(nèi)天體可能會(huì)分散在多個(gè)層級(jí)的不同星表子集文件中，此時(shí)需要遍歷該區(qū)域的父節(jié)點(diǎn)和子節(jié)點(diǎn)查詢。

　　如果想要在切分后星表中檢索到某個(gè)區(qū)域，例如第 k 層，HEALPix 編號(hào)為 pix 的節(jié)點(diǎn)所覆蓋區(qū)域內(nèi)全部天體。首先從該 pix 節(jié)點(diǎn)開始，向上遞歸訪問其父節(jié)點(diǎn)，直到訪問到第一個(gè)存在星表子集文件的節(jié)點(diǎn)，讀取文件內(nèi)所有滿足編號(hào)在目標(biāo)區(qū)域內(nèi)的天體到內(nèi)存中。若 k 小于 kₘₐₓ，則還需要依次遍歷該 pix 節(jié)點(diǎn)的所有子節(jié)點(diǎn)，并將子節(jié)點(diǎn)中存在對應(yīng)的星表子集的文件全部讀取，最終得到符合要求的天體信息。

　　若待檢索區(qū)域是一個(gè)多邊形，則可以將該多邊形區(qū)域?qū)氲?MOC 中，得到對應(yīng)區(qū)域的 UNIQ 列表。再對該 UNIQ 列表中的每個(gè)區(qū)域進(jìn)行上述的區(qū)域查詢，并匯總，即可得到最終的檢索結(jié)果。

　　3.2 對坐標(biāo)的檢索

　　若給定一個(gè)坐標(biāo)和閾值，要求檢索以坐標(biāo)為圓心，閾值為半徑內(nèi)的所有目標(biāo)，首先使用 HEALPix 自帶的 querydisc 函數(shù)計(jì)算得到在最大層次下被該區(qū)域覆蓋得到的 pix 列表，將列表中的每一個(gè) pix 執(zhí)行上述區(qū)域查詢，并對已經(jīng)讀取的星表做好標(biāo)記避免重復(fù)讀取。

　　對于兩個(gè)在天球上的坐標(biāo) (α₁,δ₁) 和 (α₂,δ₂)，他們之間的角距離 d 可由 Haversine 公式給出：

　　d=2arcsin√[sin²((δ₁-δ₂)/2)+cosδ₁cosδ₂sin²((α₁-α₂)/2)]

　　對內(nèi)存中的天體分別計(jì)算并將與坐標(biāo)的角距離小于閾值的目標(biāo)加入到返回列表中，得到檢索的結(jié)果。

　　而在最近鄰檢索時(shí)，待檢索坐標(biāo)與最近鄰天體可能不在 kₘₐₓ層級(jí)下同一個(gè)區(qū)域中，所以需要讀取與該區(qū)域相鄰接區(qū)域內(nèi)的所有天體。讀取后通過 Haversine 公式計(jì)算與待檢索坐標(biāo)的距離，找到與待檢索目標(biāo)距離最近的天體。

　　4 數(shù)據(jù)的高效存取

　　數(shù)據(jù)在內(nèi)存中通常是以對象、結(jié)構(gòu)體、數(shù)組、列表等形式連續(xù)或離散地存儲(chǔ)于不同的區(qū)域供 CPU 訪問和操作。當(dāng)內(nèi)存中的數(shù)據(jù)需要被寫入到文件或通過網(wǎng)絡(luò)傳輸時(shí)，則需要對其進(jìn)行編碼(也稱之為序列化 (serialization)、編組 (marshalling))為字節(jié)序列，反過程則稱之為解碼(或反序列化 (deserialization)、反編組 (unmarshalling))。例如將內(nèi)存中的天體以 CSV 格式保存到硬盤即是對星表序列化的過程，對 CSV 星表的讀取則是反序列化過程。

　　在切分過程中，星表的每個(gè)天體都會(huì)被讀取兩次，寫入兩次，數(shù)十億天體將會(huì)產(chǎn)生大量硬盤 I/O 并消耗大量的時(shí)間。一方面，由于星表的坐標(biāo)等數(shù)據(jù)都是以浮點(diǎn)數(shù)的形式存儲(chǔ)在內(nèi)存中，如果簡單地以 CSV 等文本格式存儲(chǔ)，需要解析每行的各個(gè)表項(xiàng)以及對浮點(diǎn)數(shù)與字符串之間相互轉(zhuǎn)換，由此產(chǎn)生的計(jì)算量和耗時(shí)也不容小覷。另一方面，在切分時(shí)，程序只讀取天體的坐標(biāo)，其他的物理參數(shù)則不必考慮。

　　為了盡可能地提升性能，本文以開源的高性能序列化組件協(xié)議緩沖區(qū) (Protocol buffers, Protobuf) 為基礎(chǔ)，并針對星表這一格式的數(shù)據(jù)定義了對應(yīng)的序列化格式，實(shí)現(xiàn)對星表數(shù)據(jù)的高效存取。

　　Protobuf 能夠?qū)?nèi)存中數(shù)據(jù)按照事先約定的格式寫入程序與程序之間，或者通過網(wǎng)絡(luò)傳遞，也可以寫入到磁盤保存。功能上類似于可擴(kuò)展標(biāo)記語言 (extensible markup language, XML) 和 JavaScript 對象簡譜 (JavaScript Object Notation, JSON)。但是與之不同的是，Protobuf 需要事先約定好一套協(xié)議，稱之為接口描述語言 (interface description language, IDL)，讀寫時(shí)需要持有相同的 IDL。數(shù)據(jù)按照 IDL 規(guī)定的格式以二進(jìn)制的形式存儲(chǔ)，并使用多種壓縮算法，避免了數(shù)字與字符串之間的相互轉(zhuǎn)換。并且支持 C++、JAVA、Python、Go 等主流的語言。

　　本文定義了星表在 Protobuf 中的傳輸格式，包括一個(gè) “Catalog” 項(xiàng)用于儲(chǔ)存星表信息，內(nèi)容如表 2 所示，其中包括星表文件的文件名、UNIQ 編號(hào)、星表頭信息、總行數(shù)信息等;同時(shí)內(nèi)部有任意數(shù)量的 “CatalogLine” 項(xiàng)，用于存儲(chǔ)天體信息，內(nèi)容如表 3 所示，其中包括了天體的坐標(biāo)、區(qū)域編號(hào)，分別以浮點(diǎn)數(shù)和整數(shù)進(jìn)行存取。而天體的其他信息則保存為字符串類型，以減少對數(shù)據(jù)的解析時(shí)間。

　　根據(jù)定義的傳輸格式，Protobuf 可以自動(dòng)生成不同編程語言下的序列化和反序列化的代碼。這使得在星表切分、歸并和檢索時(shí)，能夠根據(jù)需要，將內(nèi)存中的星表序列化為 Protobuf 格式并保存到硬盤對應(yīng)的文件中，以及從硬盤中讀取文件并反序列化回內(nèi)存，實(shí)現(xiàn)了對星表數(shù)據(jù)的高效存取。

　　5 多個(gè)大規(guī)模星表的融合

　　在暫現(xiàn)源搜尋時(shí)往往需要同時(shí)對多個(gè)不同波段的星表，或者多個(gè)不同觀測時(shí)間星表進(jìn)行檢索。為了加快檢索速度，可以提前將多個(gè)星表合并，即對一個(gè)星表中的天體與其他星表交叉證認(rèn)，找到對應(yīng)天體并對各項(xiàng)參數(shù)加以融合，從而實(shí)現(xiàn)一次檢索就可以獲得多個(gè)星表中的數(shù)據(jù)。

　　由于不同望遠(yuǎn)鏡觀測時(shí)的精度差異等多種原因?qū)е孪嗤脑丛诓煌男潜碇械奈恢貌皇峭耆嗤摹Ｒ话銇碚f，對于誤差半徑分別為 r₁、r₂的兩個(gè)星表，當(dāng)兩個(gè)天體的角距離 d<3√(r₁²+r₂²) 的時(shí)候，就可以認(rèn)為是同一天體。位于 HEALPix 劃分區(qū)域邊緣的天體誤差半徑可能會(huì)覆蓋到相鄰區(qū)域，所以在計(jì)算某一個(gè)區(qū)域的距離時(shí)還需要讀取與該區(qū)域相鄰的幾個(gè)區(qū)域內(nèi)的天體，此時(shí)可以通過選擇合適的遍歷順序和增加緩存以提升效率。

　　5.1 星表在內(nèi)存中的索引

　　將天體按切分到不同的星表文件時(shí)，天體在文件內(nèi)的順序即切分時(shí)的讀入順序。若需要檢索某個(gè)文件內(nèi)的天體時(shí)，順序檢索效率不高，可以將該文件全部讀入到內(nèi)存，并在內(nèi)存中建立 KD-Tree，實(shí)現(xiàn)二維檢索，提高檢索速度。

　　KD-Tree 是一種對空間多維數(shù)據(jù)高效索引的數(shù)據(jù)結(jié)構(gòu)。KD-Tree 會(huì)在建立時(shí)在每層選擇不同維度的數(shù)據(jù)計(jì)算中位數(shù)并作為索引節(jié)點(diǎn)，最終構(gòu)建出一顆比較均衡的二叉樹。而 KD-Tree 在檢索時(shí)，則按照層數(shù)分別將不同維度的數(shù)據(jù)與樹中的數(shù)據(jù)進(jìn)行對比，從而能夠高效檢索高維數(shù)據(jù)。

　　根據(jù)之前的切分情況，以星表文件作為建立 KD-Tree 的單位，若某次檢索的區(qū)域覆蓋了多個(gè)星表文件，則需要在內(nèi)存中建立多個(gè) KD-Tree。此時(shí)可以將建立好的 KD-Tree 緩存到內(nèi)存中，供之后的檢索使用，避免重復(fù)讀取，提高效率。而緩存大小不可能無限大，特別是在同時(shí)運(yùn)行多個(gè)任務(wù)時(shí)，內(nèi)存中的 KD-Tree 數(shù)量達(dá)到上限，則可以選擇合適的緩存替換算法，例如先進(jìn)先出 (FIFO)，或者最近最久未使用 (LRU) 等，替換緩存中的 KD-Tree。各種替換算法在計(jì)算機(jī)等領(lǐng)域已經(jīng)被深入研究過，本文將通過實(shí)驗(yàn)探究不同緩存算法在星表證認(rèn)時(shí)的性能差異以及最優(yōu)緩存的大小。

　　5.2 使用 Peano-Hilbert 曲線遍歷星表

　　在對 HEALPix 按照 NESTED 模式編號(hào)順序遍歷時(shí)，路徑在空間表現(xiàn)為 Z 曲線的形式。該順序可以高效地為平面或空間中的一組點(diǎn)集構(gòu)建四叉樹或者八叉樹。而 Z 曲線是在不同區(qū)域的邊界處有比較大的轉(zhuǎn)折，即會(huì)跳轉(zhuǎn)到一個(gè)不相鄰的區(qū)域中。在對星表遍歷時(shí)往往會(huì)導(dǎo)致之前的緩存大量失效，降低效率。為了克服這個(gè)缺點(diǎn)，本文嘗試在對星表遍歷時(shí)使用 Peano-Hilbert 曲線代替原有的 Z 曲線。該曲線通過對 HEALPix 的區(qū)域重新進(jìn)行編號(hào)，實(shí)現(xiàn)了相鄰編號(hào)的兩個(gè)區(qū)域在天球上也相鄰，因此對星表進(jìn)行遍歷時(shí)會(huì)一直訪問相鄰的區(qū)域，能夠充分利用已有的緩存。按編號(hào)順序分別對 NESTED 自帶編號(hào)和 Peano-Hilbert 編號(hào)對星表進(jìn)行遍歷。

　　6 實(shí)驗(yàn)分析

　　本實(shí)驗(yàn)使用 LAMOST DR7 星表作為測試數(shù)據(jù)，運(yùn)行平臺(tái)使用國家天文科學(xué)數(shù)據(jù)中心提供的云虛擬機(jī)，配置如下：

　　CPU 為 Intel (R) Xeon (R) CPU E5-2690 4 Core @ 2294.686 MHz;

　　內(nèi)存為 8 GB;

　　硬盤 I/O 速度為 82.5 MB/s;

　　操作系統(tǒng)為 Ubuntu 20.04.4 LTS;

　　編程語言為 C++、Python。

　　6.1 不同切分方式產(chǎn)生的星表文件對比

　　分別對測試星表進(jìn)行 “激進(jìn)”、“中立”、“保守” 模式的動(dòng)態(tài)切分和按固定層級(jí)的切分。測試動(dòng)態(tài)切分使用的切分閾值為 1000。

　　可以看出，“激進(jìn)” 方式切分出的星表的天體數(shù)量均小于閾值，大部分集中在 [200,1000] 之間。“中立” 方式切分出的星表的天體數(shù)量絕大部分集中在 [1000,4000] 之間;由于某些南天區(qū)的目標(biāo)過少，導(dǎo)致極少量的文件小于閾值;同時(shí)，“中立” 模式下切分出的星表的變異系數(shù)也是所有模式中最低的，說明星表的離散程度小于其他模式。“保守” 方式切分出的星表均大于 1000，特別是一部分星表的天體數(shù)量過大，最大約 10⁶，其變異系數(shù)也是各個(gè)方式中最大的。三個(gè)固定層級(jí)切分的結(jié)果變異系數(shù)相差不大，皆大于 “中立” 和 “激進(jìn)” 模式，小于 “保守模式｡

　　切分出的結(jié)果與預(yù)期相符，“中立” 方式切分出的星表較其他兩種切分方式和固定切分更為均勻｡

　　6.2 Protobuf 格式與格式星表的讀取速度對比

　　本文將 LAMOST DR7 星表分別在 “中立” 模式下按照不同的閾值切分為 CSV 格式和 Protobuf 格式，并分別對其全部讀取。可以看出，隨著切分閾值的增加，兩種格式的星表讀取時(shí)間都逐漸減小直至達(dá)到穩(wěn)定｡Protobuf 格式在文件體積略大于 CSV 格式的情況下，當(dāng)閾值大于 1000 之后的星表讀取速度比 CSV 格式快 50% 以上，基本達(dá)到硬盤連續(xù)讀取 I/O 的速度｡之后，本文又對這兩種格式的星表進(jìn)行了自交叉證認(rèn) (self-matching) 的測試，即先將未切分星表全部讀入到內(nèi)存并按 HEALPix 編號(hào)排序作為源星表，按天體在源星表中的順序計(jì)算對應(yīng)切分后目標(biāo)星表文件的編號(hào)｡若文件已存在于內(nèi)存中，則與內(nèi)存中的 KD-Tree 交叉證認(rèn);若未存在，則刪除內(nèi)存中的 KD-Tree，并讀取對應(yīng)的文件，建立新的 KD-Tree 進(jìn)行交叉證認(rèn)｡統(tǒng)計(jì)不同閾值下的證認(rèn)時(shí)間｡可以看出，隨著閾值的增大，兩種格式的星表證認(rèn)時(shí)間都先減小后增大｡CSV 格式的星表在閾值為 1000 左右時(shí)速度最快，需要 155s 左右;而 Protobuf 格式的星表在閾值為 10000 時(shí)速度最快，只需約 116s 即可完成｡因此 Protobuf 格式對星表存儲(chǔ)和交叉證認(rèn)時(shí)的效率優(yōu)于 CSV｡

　　6.3 不同曲線與緩存算法在星表證認(rèn)時(shí)的效率對比

　　對不同切分層級(jí)，分別用 HEALPix NESTED 模式下自帶的 Z 曲線和 Peano-Hilbert 曲線遍歷整個(gè)天區(qū)，在遍歷到第zˉ個(gè)區(qū)域時(shí)，同時(shí)讀取z˙和以zˉ為中心的上、下、左、右、左上、左下、右上、右下共九個(gè)區(qū)域的信息｡若內(nèi)存中已經(jīng)緩存了該區(qū)域的信息，則為緩存命中;若不存在，則為緩存失效｡此時(shí)將該區(qū)域的信息從磁盤中讀取并加入到緩存，若緩存區(qū)已滿，則調(diào)用替換算法選擇一個(gè)區(qū)域替換，分別使用最久未使用 (LRU) 替換策略與先入先出 (FIFO) 緩存替換策略進(jìn)行測試｡

　　記錄在不同層級(jí)和不同緩存大小的緩存失效次數(shù)，兩種曲線和兩種緩存算法的四種組合在四個(gè)不同的層級(jí)下的結(jié)果｡可以看出：緩存越大，緩存失效的次數(shù)會(huì)越小，但單位緩存增加帶來的性能提升也隨之減少｡LRU 的效率要好于 FIFO，Peano-Hilbert 曲線的效率要好于 Z 曲線｡遍歷曲線和緩存替換算法的效率差距幾乎不受 k 值大小的影響｡

　　少量的緩存大小就可以大幅度減少緩存失效的次數(shù)，提高效率｡例如使用 Peano-Hilbert+LRU 方法，與不緩存 + FIFO 相比，當(dāng)緩存大小為 9 時(shí)，就可以減少大約 66.7% 的讀取次數(shù);緩存大小增大到 20 時(shí)，則可以減少 77.4% 的讀取次數(shù)｡

　　聯(lián)合使用 Peano-Hilbert 曲線和 LRU 算法與 Z 曲線和 FIFO 之間的差異在緩存大小為 9 時(shí)達(dá)到最大，接近 1 倍｡

　　由此可見，通過使用 Peano-Hilbert 曲線和 LRU 的緩存替換策略可以有效提升在星表融合時(shí)對天體交叉證認(rèn)的效率｡

　　7 結(jié)論

　　針對大規(guī)模星表的高效檢索問題，本文以文件作為存儲(chǔ)的基礎(chǔ)，提出了包括星表的切分、索引、合并在內(nèi)的一系列解決方案｡一種自適應(yīng)密度的切分算法，能夠在有限內(nèi)存的情況下將大規(guī)模星表切分成多個(gè)限定大小范圍的子星表，方便后續(xù)的檢索等操作｡經(jīng)過測試，本文提出的 “中立” 模式切分算法可以將切分后的星表大小限制在[T,4×T]之間，比其他模式更均勻，為星表高效檢索奠定了基礎(chǔ)｡

　　為了更快地檢索星表，本文引入開源序列化庫 Protocol Buffers 作為星表在磁盤中的存儲(chǔ)格式｡該庫能夠?qū)⑻祗w坐標(biāo)以二進(jìn)制格式存儲(chǔ)，在切分與檢索時(shí)避免了浮點(diǎn)數(shù)與字符串之間的相互轉(zhuǎn)換造成的性能消耗，提升了星表的讀取效率｡

　　在合并星表時(shí)，需要讀取相鄰星表以解決漏源問題｡本文通過改進(jìn) HEALPix 的遍歷方式，以 Peano-hilbert 曲線代替原有的 Z 曲線，并引入緩存和替換算法，減少了星表合并時(shí)的文件讀取次數(shù)，提高了合并的速度｡

　　在本文工作的基礎(chǔ)上，未來可對諸多部分進(jìn)行進(jìn)一步優(yōu)化與創(chuàng)新｡例如，本文所述檢索方法與測試結(jié)果是在單線程下取得，未來可將該方法擴(kuò)展至多線程，或利用分布式集群等方式將檢索任務(wù)發(fā)送到多個(gè)子節(jié)點(diǎn)處理｡另外，本文交叉證認(rèn)的部分僅依靠坐標(biāo)間角距離來判斷是否同源，未來可將天體的亮度或其他參數(shù)作為判據(jù)引入，以提升在密集星場等復(fù)雜場合交叉證認(rèn)的準(zhǔn)確度｡本文針對文件而設(shè)計(jì)諸多創(chuàng)新技術(shù)，未來可以嘗試移植到數(shù)據(jù)庫等其他非文件的場合，并加以優(yōu)化，以期待實(shí)現(xiàn)更好的性能｡

張琦乾;樊東衛(wèi);崔辰州,中國科學(xué)院國家天文臺(tái);中國科學(xué)院大學(xué);國家天文科學(xué)數(shù)據(jù)中心,202303