楊樂強，王建立，姚凱男，李宏壯，陳 璐，邵 蒙

（中國科學院 長春光學精密機械與物理研究所，吉林 長春 130033）

1 引 言

地基大口徑光學望遠鏡的成像分辨率受到大氣湍流引起的波前畸變影響會下降。自適應光學系統(tǒng)通過探測并補償入射波前的相位畸變，能夠有效克服大氣湍流引起的相位畸變影響，目前已成為提高地基大口徑光學望遠鏡成像分辨率的必要手段［1-3］。在自適應光學系統(tǒng)中，波前處理器負責解算入射波前的畸變信息，根據(jù)控制算法計算產(chǎn)生相應的驅動控制信號，控制波前校正器補償入射波前畸變，是自適應光學系統(tǒng)的運算核心，其處理能力與計算延時直接決定了自適應光學系統(tǒng)的閉環(huán)動態(tài)性能。

隨著地基望遠鏡口徑的逐步增大，自適應光學系統(tǒng)波前校正單元數(shù)也不斷增加。為滿足4 m級光學望遠鏡的波前探測需求，自適應光學系統(tǒng)的規(guī)模需要提高至千單元量級，相應的波前處理器的計算規(guī)模隨之呈平方倍增加。另一方面，波前探測器的采樣頻率要超過1 000 Hz，才能保證自適應光學系統(tǒng)具有足夠的誤差抑制帶寬，從而抑制大氣湍流引起的動態(tài)波前畸變。因此，大規(guī)模自適應光學系統(tǒng)對波前處理器的數(shù)據(jù)吞吐規(guī)模、計算速度提出了很高的要求［4-6］。

為了降低系統(tǒng)延遲，提高數(shù)據(jù)處理速度，傳統(tǒng)的自適應光學系統(tǒng)大多采用基于數(shù)字信號處理器（Digital Signal Processor，DSP）陣列構成的波前處理器以及基于現(xiàn)場可編程門陣列（Field Programmable Gate Array，F(xiàn)PGA）的波前處理器［7-8］。美國星火靶場的SOR3.5 m 望遠鏡自適應光學系統(tǒng)，系統(tǒng)規(guī)模為941 單元，采用9 塊FPGA 級聯(lián)作為波前處理器，計算延時約為0.297 ms［9-10］。2012 年，Veran 針對三十米望遠鏡（Thirty Meter Telescope，TMT）在近紅外波段觀測的多層共軛自適應光學系統(tǒng)的計算需求，提出了采用16 塊FPGA 刀片的ATCA 架構的處理方案［11］。中國科學院成都光電所的61 單元自適應光學系統(tǒng)，其波前處理器采用17 個多核DSP 組成，峰值計算速度可達8.5 億次/秒，實驗結果表明，該波前處理器的運算延時為330 μs，誤差抑制帶寬可達54 Hz［12］。長春光機所研制的基于FPGA 的千單元波前處理板卡，采用10 塊板卡組成波前處理器，采用多線并行流水算法縮短波前處理延時，提高系統(tǒng)控制帶寬［4-5］。2016 年，成都光電所報道了1 m 太陽望遠鏡自適應光學系統(tǒng)，該系統(tǒng)采用151 單元變形鏡，波前處理器采用FPGA 與DSP 協(xié)同工作的方式，能夠在系統(tǒng)采樣頻率為3 100 Hz 的條件下實現(xiàn)穩(wěn)定閉環(huán)，誤差抑制帶寬可達170 Hz［13-14］。

基于DSP 及FPGA 的波前處理器，具有并行程度高、計算數(shù)據(jù)快以及集成度高等優(yōu)點，但其多芯片分布式的計算架構調試困難，擴展性較差。隨著商業(yè)化圖形處理器（Graphics Processing Unit，GPU）與多核中央處理器（Central Processing Unit，CPU）的并行計算能力的逐步增強，自2010 年 以 來，GPU 與 多 核CPU 逐 漸 成 為大口徑望遠鏡自適應光學系統(tǒng)波前處理計算架構的研究熱點［15-17］。2012 年，Basden 等針對歐洲南方天文臺計劃建設的極大型望遠鏡系統(tǒng)中的自適應光學系統(tǒng)波前處理要求，采用CPU 及GPU 方案對其波前處理能力進行了分析［18］。對于規(guī)模為40×40 的自適應光學系統(tǒng)，采用E5645的雙CPU 系統(tǒng)的波前處理延時約為323 μs。Basden 等還分析了基于GPU 的波前處理架構能力及其適應性，并提出了一系列優(yōu)化措施。2015年，Basden 等針對規(guī)模數(shù)為80×80 的自適應光學系統(tǒng)，采用基于GPU 的波前處理方案，實現(xiàn)了550 Hz 的波前處理速度［19-20］。2018 年，該課題組采用多核多節(jié)點CPU 處理器作為波前處理器處理同樣規(guī)模的自適應光學系統(tǒng)，整個鏈路延時為1.35 ms，處理延時測量結果的方差為20 μs［21-22］。

上述GPU 波前處理架構研究大多基于對計算架構的波前處理計算速度及其適應性進行分析，僅停留在分析和模擬的階段，實際的大規(guī)模自適應光學系統(tǒng)的系統(tǒng)級性能分析結果相對較少。本文在基于GPU 的波前處理架構的基礎上，通過桌面961 單元自適應光學系統(tǒng)的動態(tài)校正實驗，驗證了基于GPU 的波前處理架構的動態(tài)像差抑制能力，通過湍流屏動態(tài)模擬分析了不同格林伍德頻率下自適應光學系統(tǒng)的成像結果，并對其動態(tài)性能進行了分析。

2 桌面961 單元自適應光學系統(tǒng)組成

自適應光學系統(tǒng)的原理示意圖如圖1 所示。該系統(tǒng)主要由波前傳感器、波前處理器以及波前校正器3 個子系統(tǒng)構成。其中，波前傳感器負責測量入射波前的畸變信息，并將測量結果發(fā)送給波前處理器；波前處理器接收波前傳感器的測量結果，并計算得到需要施加到波前校正器上的控制量；波前校正器根據(jù)控制結果產(chǎn)生相應波前補償入射的畸變波前，完成自適應光學系統(tǒng)波前控制閉環(huán)。經(jīng)過自適應光學系統(tǒng)校正后的遠場成像，可以達到近衍射極限的成像分辨率。

圖1 自適應光學系統(tǒng)原理Fig.1 Principle diagram of adaptive optics system

圖2 為桌面961 單元自適應光學系統(tǒng)光路示意圖。白光光源發(fā)出的光束經(jīng)過湍流模擬器后通過準直透鏡準直為平行光束，平行光束入射到兩級快速反射鏡構成的傾斜像差校正系統(tǒng)，進行傾斜像差校正，經(jīng)過快速反射鏡的光束再經(jīng)過擴束元件擴大后入射到961 單元變形鏡上進行高階像差的閉環(huán)校正。經(jīng)過變形鏡校正后的光束通過分色鏡將光分為兩束，其中波段為700～900 nm 的光經(jīng)過透鏡會聚成像到高分辨率成像相機中，而另外一束500～700 nm 的光則經(jīng)過縮束系統(tǒng)將光束口徑縮小后入射到夏克-哈特曼波前傳感器進行波前探測。夏克哈特曼波前傳感器與變形鏡構成閉環(huán)控制系統(tǒng)進行系統(tǒng)的入射波前殘差校正。系統(tǒng)中961 單元變形鏡的有效光瞳口徑為235 mm，變形鏡的促動器采用PZT 型壓電陶瓷促動器，單個促動器行程為5 μm，促動器間距為7 mm，呈矩形排布方式。夏克-哈特曼波前傳感器的有效子孔徑數(shù)為1 020，微透鏡陣列排布為37×37，單個微透鏡口徑為200 μm，焦距為7 mm，微透鏡陣列成像點斑經(jīng)過中繼匹配鏡組縮放后入射到波前傳感器相機靶面中，以滿足子孔徑成像艾里斑的尺寸要求；波前傳感器相機采用FirstLight 公司Ocam2 相機，相機靶面分辨率為240×240，像元尺寸為24 μm，最高采樣幀頻可達2 000 Hz。成像相機采用濱松公司的ORCA Flash4.0 科學級成像相機，相機分辨率為2 048×2 048，像元尺寸為6.5 μm。

圖2 961 單元自適應光學系統(tǒng)測試光路Fig.2 Optical layout of 961-element adaptive optics testbed system

3 基于GPU 的波前處理器的設計與優(yōu)化

3.1 波前處理算法

自適應光學系統(tǒng)波前處理算法主要由波前斜率解算、波前復原運算以及波前控制量解算三部分構成，其計算流程如圖3 所示。波前傳感器獲取的帶有入射波前畸變信息的圖像，并將圖像傳輸給波前處理器，波前處理器接收波前傳感器的圖像首先進行波前斜率計算，夏克哈特曼波前傳感器第k個子孔徑波前斜率的計算公式如下：

其中：（xi，yj）是像素在子孔徑的x和y方向上的坐標，Ii，j則是在（xi，yj）坐標點上的像素灰度值。（xc，yc）為實際入射波前經(jīng)過夏克哈特曼波前傳感器單個微透鏡后聚焦在探測器靶面上的光斑實際質心位置，（xref，yref）則為入射波前為理想平面波時探測器靶面形成的光斑質心位置。由實際質心位置與理想光斑質心位置的偏差，就可以獲得夏克-哈特曼波前傳感器內子孔徑的二維波前斜率信息（Δxi，Δyi）。將每個子孔徑內的斜率向量按x和y方向重新排列成一維列向量，就可以獲得波前斜率向量：

基于波前斜率向量，本文通過直接斜率法進行波前復原計算。該方法通過矩陣向量乘法的方式得到相應的復原電壓，可以表示為：

上述形式的矩陣表示為：

其中：E是誤差向量，是復原運算的結果；D是響應矩陣的廣義逆矩陣，一般稱為控制矩陣。

根據(jù)波前誤差向量進行波前控制運算，本文采用自適應光學系統(tǒng)中常用的PI 控制運算［12，23］。

3.2 基于GPU 的波前處理架構

基于GPU 的波前處理架構如圖3 所示，主要由圖像采集卡、負責波前處理的GPU 以及負責協(xié)同控制與數(shù)據(jù)傳輸?shù)腃PU 構成。在波前處理過程中，哈特曼波前探測器作為數(shù)據(jù)源，將微透鏡陣列所成的光斑陣列圖像通過CameraLink 線纜發(fā)送到GPU 所在的波前處理計算機的圖像采集卡中。圖像采集卡在接收完一幀圖像信息后，通過DMA 的方式將波前探測器的圖像傳輸?shù)缴衔粰C內存當中，同時以中斷的形式通知CPU 波前探測器；CPU 響應該中斷，將接收到的圖像數(shù)據(jù)通過PCIe 總線傳輸?shù)紾PU 顯存部分。GPU接收到圖像后，首先進行波前斜率計算即子孔徑質心運算，得到質心偏移向量；然后進行波前復原運算，將得到的質心偏移向量和初始化過程中就已經(jīng)拷貝到顯存當中的波前復原矩陣相乘，得到波前誤差向量；最后根據(jù)波前誤差向量，進行波前控制運算，再根據(jù)PID 控制參數(shù)計算得到變形鏡驅動電壓控制量。GPU 同樣通過PCIe 總線傳輸?shù)姆绞綄@存中的驅動電壓控制向量傳回CPU 內存當中，CPU 將驅動電壓控制量打包并使用UDP 通信協(xié)議通過千兆網(wǎng)口將它發(fā)送到數(shù)字模擬轉換卡（DA），轉換后的模擬控制量經(jīng)過放大器放大后驅動變形鏡對應的壓電陶瓷促動器產(chǎn)生位移形變，生成共軛波前，完成自適應光學系統(tǒng)閉環(huán)校正的一次迭代過程。

圖3 基于GPU 的自適應光學系統(tǒng)波前處理算法流程Fig.3 Flow chart of wavefront processor of GPU-based adaptive optics system

3.3 基于GPU 的波前處理算法優(yōu)化

在GPU 并行計算框架中，CUDA 和OPENCL 是兩種常見的GPU 并行計算框架。由于CUDA 計算框架廣泛應用于深度學習、大規(guī)模并行計算中，本文選擇CUDA 計算框架作為GPU的并行計算框架。

為了提高基于GPU 波前處理架構的波前處理速度，本文在基于CUDA 計算框架的波前處理算法中進行了一系列的算法映射與優(yōu)化設計，盡可能地提高GPU 計算資源的利用效率，進而提高波前處理速度。主要的算法映射過程與優(yōu)化處理步驟如下：

（1）減少無效圖像帶來的傳輸開銷。在基于GPU 的波前處理機構中，圖像數(shù)據(jù)在內存與顯存之間的交換導致GPU 計算延時，需要盡可能降低由傳輸數(shù)據(jù)帶來的傳輸延時。夏克-哈特曼波前傳感器在波前斜率解算過程中，其有效數(shù)據(jù)僅為子孔徑圖像信息，而其余部分為無效圖像數(shù)據(jù)，因此可以通過去除圖像中無效像素信息的方式，降低內存與GPU 之間傳輸?shù)耐ㄐ艛?shù)據(jù)量。根據(jù)961 單元自適應光學系統(tǒng)夏克-哈特曼波前傳感器的設計參數(shù)，有效子孔徑占據(jù)的像素數(shù)為1020×6×6=36 700，而相機全靶面像素數(shù)為240×240=57 600，通過無效像素排除，可以降低36%的數(shù)據(jù)傳輸量。

（2）波前斜率計算算法映射與優(yōu)化。由于夏克-哈特曼子孔徑之間的斜率解算具有天然的并行計算特點，同時每個子孔徑僅有36 個像素參與計算，計算規(guī)模較小，因此本文將每個子孔徑的斜率計算分配到GPU 的每個計算線程當中，1 020 個子孔徑的斜率計算對應在并行的1 020個GPU 線程當中。由于計算規(guī)模較小，在質心計算過程中，GPU 的內存訪問開銷不可忽略，為了提高質心計算速度，本文將夏克-哈特曼波前傳感器的圖像按照子孔徑排序的方式進行了重排與向量化操作，使每個并行計算的線程能夠快速連續(xù)訪問子孔徑的圖像信息，提高顯存的讀取效率。在得到子孔徑質心向量后，需要對質心向量解算平均值，將傾斜像差數(shù)據(jù)單獨發(fā)送到快速反射鏡驅動器進行傾斜像差校正，減小變形鏡的行程負擔。GPU 計算得到的質心誤差向量是中間結果，位于訪問速度最快的寄存器當中，因此質心誤差的平均值計算采用CUDA 架構中的__shuffle__指令，通過該指令線程間可以進行信息交互，提高內存訪問效率進而提高計算速度。

（3）波前復原運算算法映射。波前復原向量通過質心誤差向量與復原矩陣相乘得到，其運算實質是向量矩陣乘法，因此本文使用CUDA 架構中的cublas 庫，采用列主序的數(shù)據(jù)排列方式進行矩陣向量乘法運算。

（4）波前控制算法映射。得到波前復原誤差向量后，波前控制運算可以在GPU 中映射為向量的乘加運算，同樣采用cublas 庫進行基于PID的控制計算得到波前控制向量。為了減少CPU與GPU 之間的通信開銷，波前控制向量的存儲采用CUDA 架構中的零拷貝內存，CPU 可以直接訪問該內存，減少了通信開銷。

4 實驗結果與分析

根據(jù)上述波前處理器架構與算法設計，將基于GPU 的波前處理器應用到961 單元變形鏡自適應光學系統(tǒng)中。其中，波前處理器采用的GPU 為GeForce RTX2080Ti 顯卡，CPU 采用Intel i9-9900，8 核心16 線程，主頻為3.6 GHz。

4.1 系統(tǒng)計算時延

為了驗證GPU 實現(xiàn)波前處理的有效性，本文還進行了基于CPU 的波前處理實驗，與GPU波前處理延時結果進行對比，其中CPU 波前處理中波前復原運算中最耗時的矩陣乘法部分采用OpenCV 庫函數(shù)進行矩陣乘法，統(tǒng)計10 000 幀波前處理延時后，得到波前處理延時均值，如表1所示。由表1 可以看出，基于CPU 的波前處理延時平均值為1.53 ms，當系統(tǒng)采樣頻率超過1 000 Hz時，無法滿足波前處理器的實時處理要求；而基于GPU的波前處理延時平均值為0.3 ms，滿足系統(tǒng)采樣頻率為2 000 Hz時的波前處理速度要求。

表1 波前處理延時結果Tab.1 Results of wavefront process delay

4.2 靜態(tài)像差校正實驗

首先對桌面自適應光學系統(tǒng)進行了靜態(tài)像差展平實驗。圖4 為一組變形鏡靜態(tài)像差校正實驗校正前后夏克-哈特曼波前傳感器測量得到的波前像差結果，由圖可知，校正前系統(tǒng)像差峰峰值（Peak Valley，PV）為5.30λ（λ=600 nm），均方根值（Root Mean Square，RMS）為0.81λ；校正后，波 前 殘 差PV 值 為0.24λ，RMS 值 為0.05λ。校正前后遠場光斑圖像如圖5 所示，可以看出，校正前遠場光斑圖像的灰度峰值為3 347，校正后的灰度峰值為45 876，遠場光斑的能量集中度明顯提高。

圖4 自適應光學系統(tǒng)校正前后系統(tǒng)波前像差Fig.4 System wavefront errors before and after optical adaptive corrections

圖5 校正前后光源遠場光斑圖像Fig.5 Far-field images of light source before and after correction

4.3 動態(tài)像差校正實驗

動態(tài)像差校正實驗通過湍流模擬器模擬大氣擾動。湍流模擬器通過轉動刻蝕在位相屏表面特定的光程差部分形成動態(tài)像差，動態(tài)像差能夠模擬大氣湍流的變化情況。其中，湍流模擬器所模擬的大氣格林伍德頻率與風速之間的關系為：

其中：r0代表大氣相干長度，vw代表風速，通過改變湍流模擬器轉速即可改變模擬的大氣格林伍德頻率。實驗中，湍流模擬器等效模擬的大氣相干長度為11 cm，湍流屏轉速為1 rad/s，等效模擬的格林伍德頻率為60 Hz。

為驗證變形鏡的動態(tài)性能，實驗中采用變形鏡進行像差校正，快速反射鏡不進行傾斜校正，僅作為普通反射鏡使用（實際應用中通常采用快速反射鏡校正行程較大的傾斜像差）。由湍流模擬器產(chǎn)生的傾斜像差全部由變形鏡進行校正，根據(jù)校正前后傾斜像差與波前RMS 值的動態(tài)變化，衡量波前處理系統(tǒng)的動態(tài)性能。首先，測試動態(tài)條件下自適應光學系統(tǒng)對傾斜像差的校正結果，如圖6 和圖7 所示。可以看出，自適應光學系統(tǒng)在校正前，由湍流屏引入的系統(tǒng)X軸傾斜量均方根值為0.18″，Y軸傾斜量均方根值為0.26″；經(jīng)過變形鏡校正后，X軸傾斜量均方根值下降到0.02″，Y軸傾斜量均方根值下降為0.03″。校正前后系統(tǒng)入射波前（去除傾斜像差）的均方根值變化如圖8 所示，校正前系統(tǒng)入射波前為0.5λ～1.3λ（λ=600 nm），均值為0.83λ；校正后系統(tǒng)波前殘差明顯降低，波前殘差均值為0.16λ，約為100 nm。

圖6 校正前后X 軸傾斜誤差曲線Fig.6 Track errors of X-axis before and after correction

圖7 校正前后Y 軸傾斜誤差曲線Fig.7 Track errors of Y-axis before and after correction

圖8 校正前后系統(tǒng)波前殘差曲線Fig.8 Residual wavefront error before and after correction

圖9～圖11 分別為X軸傾斜、Y軸傾斜以及高階像差均方根值校正前后的功率譜變化情況。由圖可以看出，系統(tǒng)的0 dB 誤差抑制帶寬為100 Hz，高頻部分存在校正后功率譜高于校正前的情況，這是因為在工程實際過程中，為了保證系統(tǒng)具有較高的動態(tài)特性，控制參數(shù)選取較大，相位裕度較小，存在放大高頻噪聲的情況［12，23］。

圖9 校正前后X 軸傾斜誤差功率譜Fig.9 Power spectra density of X-tilt before and after correction

圖11 校正前后系統(tǒng)波前殘差功率譜Fig.11 Power spectra density of residual wavefront error before and after correction

圖10 校正前后Y 軸傾斜誤差功率譜Fig.10 Power spectra density of Y-tilt before and after correction

圖12 是不同格林伍德頻率下，961 單元自適應光學系統(tǒng)校正前后成像相機遠場能量的三維分布?？梢钥闯?，當自適應光學系統(tǒng)閉環(huán)時，遠場能量峰值得到很大提高，但隨著格林伍德頻率的提高，自適應光學系統(tǒng)校正后的效果逐漸變差，其原因是波前處理器的閉環(huán)控制帶寬有限，隨著格林伍德頻率的提高，時域誤差也隨之提高，并逐漸成為波前殘余誤差的主要成分，從而降低自適應光學系統(tǒng)的閉環(huán)控制性能。

圖12 不同格林伍德頻率下校正后遠場光斑能量分布Fig.12 Far-field images of light source before and after correction at different Greenwood frequencies

5 結 論

本文根據(jù)961 單元自適應光學系統(tǒng)波前處理器的性能要求，采用NVIDIA 公司GeForce RTX2080Ti 顯卡作為波前處理運算核心，實現(xiàn)了961 單元自適應光學系統(tǒng)的閉環(huán)波前處理運算。實驗結果表明，對于961 單元自適應光學系統(tǒng)，基于GPU 的波前處理器在系統(tǒng)采樣頻率為1 500 Hz 時，系統(tǒng)誤差抑制帶寬可達100 Hz，滿足千單元級自適應光學系統(tǒng)的實時波前處理需求。