利用反卷積方法抽取多光纖光譜的一維譜的制作方法
【專利摘要】本發(fā)明公開了一種根據(jù)二維多光纖光譜圖像形成原理，用反卷積方法來抽取可靠的一維光譜的方法。包括光纖光譜中心跡線的獲得，每條光纖光譜上所有PSF的獲得，目標(biāo)函數(shù)的構(gòu)造，噪聲的抑制方法和目標(biāo)函數(shù)的串行和并行求解算法。利用強(qiáng)的連續(xù)光譜獲取光譜中心跡線，利用單根發(fā)射線光譜來獲取所有波長(zhǎng)處的PSF；根據(jù)噪聲和光譜情況構(gòu)造目標(biāo)函數(shù)；給出了目標(biāo)函數(shù)的串行和并行求解方法。本發(fā)明能糾正經(jīng)典抽譜方法的失真，也解決反卷積方法由于存儲(chǔ)量巨大而無法計(jì)算的問題，能得到更為真實(shí)的一維譜。同時(shí)該計(jì)算方法很靈活，內(nèi)存用量少，擴(kuò)展性強(qiáng)，容易用于并行計(jì)算，也是一個(gè)計(jì)算稀疏矩陣的新算法。
【專利說明】
利用反卷積方法抽取多光纖光譜的一維譜
技術(shù)領(lǐng)域
[0001] 本發(fā)明涉及多光纖光譜的一維譜抽取技術(shù)領(lǐng)域，和大規(guī)模稀疏矩陣串行和并行計(jì) 算領(lǐng)域。
【背景技術(shù)】
[0002] 在天文領(lǐng)域中，為了同時(shí)獲取盡可能多的天體目標(biāo)的光譜，都采用了多光纖光譜 技術(shù)。也就是，每一條光纖對(duì)準(zhǔn)一個(gè)天體目標(biāo)，多個(gè)天體的光被各自的光纖引導(dǎo)到光譜儀， 然后在CCD上色散成二維光譜。采用多光纖光譜技術(shù)，一張 CCD可以同時(shí)記錄幾百條光譜。 這些光譜是以二維譜圖像的形式被CCD記錄下來的，而用于數(shù)據(jù)分析的光譜一般是一維 譜，即橫坐標(biāo)是波長(zhǎng)，縱坐標(biāo)是流量。從而需要一個(gè)算法，把CCD記錄的二維譜圖像，轉(zhuǎn)化成 可以用于分析的一維譜。
[0003] 目前流行的方法是把垂直于每一條二維光譜色散方向的每一行的流量按照一定 加權(quán)規(guī)則加起來，從而獲得一維譜。這就要求加權(quán)規(guī)則在每一行甚至所有行都是固定不變 的，但這個(gè)假設(shè)是無法成立的。
[0004] 從另一方面看，C⑶記錄的二維譜的圖像，是進(jìn)入光學(xué)系統(tǒng)之前的一維譜卷積上光 學(xué)系統(tǒng)的點(diǎn)源擴(kuò)展函數(shù)（PSF)加上各種噪聲得到的。若在二維譜圖像上，每條光纖光譜上 所有位置處的PSF已知，則可以用二維譜圖像反卷積掉已知的PSF來獲得的一維譜。這種 抽譜方法更符合物理過程，從而更為可靠。

【發(fā)明內(nèi)容】

[0005] 本發(fā)明的目的在于根據(jù)二維多光纖光譜圖像形成原理，用反卷積方法來抽取更為 可靠的一維光譜，從而矯正當(dāng)前抽譜方法的失真。
[0006] 為實(shí)現(xiàn)上述目的，本發(fā)明提出一種用反卷積方法從光纖光譜的二維CCD圖像中抽 取一維光譜的方法。該方法包括：光纖光譜中心跡線的獲得，每條光纖光譜所有波長(zhǎng)處PSF 的獲得，目標(biāo)函數(shù)的構(gòu)造，噪聲的抑制方法、目標(biāo)函數(shù)的串行和并行求解算法。
[0007] 所述的二維光譜中心跡線的獲得方法是，從較強(qiáng)的平場(chǎng)定標(biāo)燈的二維光譜，利用 重心法，求出每一條光纖光譜的每一行的重心位置，然后把每一行重心用低階多項(xiàng)式擬合， 得到每一條光纖光譜的中心軌跡。
[0008] 所述的PSF獲取方法是，對(duì)于波長(zhǎng)定標(biāo)燈的每一條光纖光譜圖像，首先獲取單根 發(fā)射線的離散輪廓。對(duì)于發(fā)射線比較弱的情況，可以延長(zhǎng)曝光來提高該發(fā)射線的信噪比。然 后對(duì)離散的輪廓進(jìn)行歸一化，就獲得了該發(fā)射線在它所處CCD位置處的離散的PSF。另外， 也可以通過激光梳來獲得，只不過這種方法比較昂貴。如果想使用光滑的PSF，可以通過利 用B樣條曲面插值離散的PSF來獲得。我們把通過這種方法獲得的PSF稱為基本的PSF。
[0009] 所述的目標(biāo)函數(shù)的構(gòu)造方法是，根據(jù)圖像的噪聲類型來構(gòu)造目標(biāo)函數(shù)。
[0010] 高斯噪聲情況：
[0011]
[0012]
[0013]
[0014] 其中，Nx表示圖像總列數(shù)；Ny表示圖像總行數(shù)；N f表示圖像記錄的光纖光譜總數(shù)； F(k，m)表示圖像上位于第k行m列的CCD像素的計(jì)數(shù)；PSF1, Jk，m)表示位于第i條光纖 光譜在第j行處的PSF在CCD上第k行第m列像素上的計(jì)數(shù)；Cl, j表示第i條光纖光譜在 第.j行處的流量值，該值是要求解的值。
[0015] 高斯噪聲下，目標(biāo)函數(shù)的矩陣形式寫法：
[0016]
C
[0017] 泊松噪聲下，目標(biāo)函數(shù)的矩陣形式寫法：
[0018]
[0019] 其中，F(xiàn)是NxXNy行的列向量，F(xiàn)(k，m)是它的第kXN x+m個(gè)元素；A是ΝχΧΝ#， Nf X Ny列的矩陣，PSF L j (k，m)是它的第k X Nx+m行，第i X Ny+j列元素；C是一個(gè)Nf X Ny行的 列向量，Cy是它第iXNy+j個(gè)元素；W是一個(gè)NxXN yRNxXNy列的對(duì)角矩陣，它的第kXNx+m 行，kXNx+m列元素是F(k，m)，其它元素為0。
[0020] 所述的噪聲抑制方法，是利用Tikhonov正則項(xiàng)來抑制噪聲?？梢允褂盟笞兞康?0、1、2次等導(dǎo)數(shù)，及這些導(dǎo)數(shù)的組合來抑制噪聲。矩陣形式的目標(biāo)函數(shù)如下：
[0021] 高斯噪聲下，目標(biāo)函數(shù)為：
[0022]
[0023]
[0024]
[0025] 其中，Γ是Tikhonov矩陣，源自導(dǎo)數(shù)的離散化；α是權(quán)重對(duì)角矩陣；Tikhonov項(xiàng) (rC) Ta (rc)對(duì)噪聲有抑制作用。
[0026] 上面的公式的推導(dǎo)，都假設(shè)噪聲是獨(dú)立同分布的。若噪聲情況復(fù)雜，則刻畫噪聲的 對(duì)稱正定矩陣W會(huì)更復(fù)雜。無論何種情況，目標(biāo)函數(shù)可以寫為如下統(tǒng)一形式：
[0027]
⑴
[0028] 所述的目標(biāo)函數(shù)的串行求解方法，是利用矩陣分塊和迭代的方式求解目標(biāo)函數(shù) (1)。由于都是稀疏矩陣，矩陣在分塊計(jì)算過程中，利用了其系數(shù)特性進(jìn)行存儲(chǔ)和計(jì)算。具 體的計(jì)算步驟如下：
[0029] 步驟0 :給定ε >〇,子矩陣塊尺度k，和向量C的長(zhǎng)度N。0 < k << N。
[0030] 步驟1 :設(shè)初始值C。: = 0, C中將要計(jì)算的變量塊的初始位置s : = 0。對(duì)于目標(biāo) 函數(shù)（1)，設(shè) R : = F，res : = RTWR，RES := res，fy = 0〇
RES-ri，轉(zhuǎn)到步驟3。
[0033] 步驟4 :若f\< ε，停止；否則f 1: = 〇,轉(zhuǎn)到步驟2。
[0034] 所述的目標(biāo)函數(shù)并行求解方法是，在計(jì)算過程中，把要求解的變量向量分成多個(gè) 變量塊，每次并行計(jì)算那些沒有交叉污染的變量塊，直到算法收斂到最小點(diǎn)。具體步驟如 下：
[0035] 步驟0 :給定ε >〇,子矩陣塊尺度k，和向量C的長(zhǎng)度N。0 < k << N。
[0036] 步驟 1 :設(shè)全局變量的初始值 C。: = 0, R : = F，res : = RTWR，RES := res，f1: = 0〇 把變量c中按照變量編號(hào)順序分為pv/^塊。除了最后一塊，每塊有k個(gè)變量。記這些塊 為= U，...，「7V/^，令p :=丨巧丨是這些矩陣塊的集合。
[0037] 步驟2 :在集合P中選取變量塊，這些變量塊在矩陣A中對(duì)應(yīng)的子矩陣的非0行號(hào) 兩兩之間都不相同。記這些變量塊的集合為A。
[0038] 步驟3 :對(duì)于每個(gè)^ e Δ，用一個(gè)計(jì)算核計(jì)算。具體如下：記:^和方分別是矩陣A和 Γ中與G相對(duì)應(yīng)列組成的子矩陣塊。令
，
[0039] 步驟4少:=W Δ，若P ,為空集，轉(zhuǎn)到步驟6 ;否則，轉(zhuǎn)到步驟2。
[0040] 步驟5:重新把變量C進(jìn)行劃分成新的pV/&l變量塊，使得每一個(gè)具有k個(gè)元 素塊與上次塊劃分的所有具有k個(gè)的元素塊至少有個(gè)變量不同。記這些塊為 = ，令0:=·^}是這些矩陣塊的集合。轉(zhuǎn)到步驟2。
[0041] 步驟6 :若f\< ε，停止；否則f 1:= 〇,轉(zhuǎn)到步驟5。
[0042] 本發(fā)明能夠利用任意形狀的PSF反卷積出任意大小的CXD光纖光譜圖像。目標(biāo)函 數(shù)的求解方法既可以用并行計(jì)算和也可以用串行計(jì)算方法進(jìn)行計(jì)算。由于符合二維光纖光 譜形成的物理過程，本發(fā)明抽取的一維譜相對(duì)于傳統(tǒng)的抽譜方法，其信噪比和分辨率會(huì)大 大提高，同時(shí)更加真實(shí)。本發(fā)明可以很自然地處理兩根相鄰的具有較嚴(yán)重相互污染的光纖 光譜圖像，這要比傳統(tǒng)方法更優(yōu)越。本發(fā)明擴(kuò)展性很強(qiáng)，可以通過設(shè)置矩陣W來刻畫不同的 噪聲，也可以通過修改Tikhonov矩陣α來根據(jù)實(shí)際情況更加靈活的對(duì)噪聲進(jìn)行抑制。
【附圖說明】
[0043] 圖1是本發(fā)明利用反卷積方法抽取多光纖光譜的一維譜的流程圖。
[0044] 圖2是本發(fā)明利用反卷積方法抽取多光纖光譜的一維譜的串行計(jì)算示意圖。
[0045] 圖3是本發(fā)明利用反卷積方法抽取多光纖光譜的一維譜的并行計(jì)算示意圖。
【具體實(shí)施方式】
[0046] 為使本發(fā)明的目的、技術(shù)方案和優(yōu)點(diǎn)更加清楚明白，以下結(jié)合具體實(shí)施例，并參照 附圖，對(duì)本發(fā)明進(jìn)一步詳細(xì)說明。
[0047] 圖1是本發(fā)明可利用反卷積方法抽取多光纖光譜的一維譜的流程圖。具體包括： 光纖光譜中心跡線的獲得，每條光纖光譜所有波長(zhǎng)處PSF的獲得，目標(biāo)函數(shù)的構(gòu)造，噪聲的 抑制方法、目標(biāo)函數(shù)的串行和并行求解算法。
[0048] 首先，我們需要知道每條光纖光譜的中心跡線，從而知道每條光纖光譜在CCD上 的位置。為此，需要有一個(gè)能發(fā)出較強(qiáng)連續(xù)譜的平場(chǎng)燈，例如碘鎢燈，照射所有光纖的輸入 端。燈光會(huì)通過光學(xué)系統(tǒng)在CCD上留下每條光纖的二維光譜圖像。利用重心法，求出每一 條光纖光譜的每一行的重心位置，然后把每一行重心用低階多項(xiàng)式擬合，得到每一條光纖 光譜的中心軌跡。
[0049] 其次，需要獲取每條光纖光譜所有波長(zhǎng)處的PSF。為此，需要有一個(gè)能發(fā)出多條單 根發(fā)射線的、且這些單根發(fā)射線盡可能均勻分布在整條光纖光譜圖像中的燈。例如，汞鎘 燈。用這種燈照射光纖輸入端，光纖會(huì)通過光學(xué)系統(tǒng)在C⑶上留下每條光纖的燈譜圖像，然 后我們可以沿著光譜中心跡線方向，根據(jù)一定的尺寸，依次提取每條光纖光譜上的所有單 根發(fā)射線。對(duì)于比較弱的單根合發(fā)射線，可以通過延長(zhǎng)曝光時(shí)間來提高該發(fā)射線的信噪比。 然后我們對(duì)離散的發(fā)射線輪廓進(jìn)行歸一化，就獲得了該發(fā)射線在它所處CCD位置處的離散 的PSF。另外，我們可以通過激光梳來獲得，只不過這種方法比較昂貴。如果想使用光滑的 PSF，可以通過利用B樣條曲面插值離散的PSF來獲得把光滑的PSF。我們把通過這種方法 獲得的PSF稱為基本的PSF，這些基本的PSF在計(jì)算過程中存儲(chǔ)在電腦內(nèi)存中。對(duì)于每一條 光纖光譜圖像，在沒有燈譜發(fā)射線的位置，可以通過與它相鄰的兩個(gè)基本PSF線性插值獲 得。這些非基本的PSF在需要的時(shí)候才計(jì)算，不存儲(chǔ)在電腦內(nèi)存中。
[0050] 然后，根據(jù)噪聲模型設(shè)定目標(biāo)函數(shù)中的W，根據(jù)光纖光譜的亮度和光纖光譜不同波 長(zhǎng)處流量相對(duì)強(qiáng)度變化，設(shè)定目標(biāo)函數(shù)中的α，從而確定了目標(biāo)函數(shù)。
[0051] 再次，需要求解目標(biāo)函數(shù)。圖2是串行求解流程圖，具體如下：
[0052] 步驟0 :給定ε > 〇,子矩陣塊尺度k，和向量C的長(zhǎng)度Ν。0 < k << Ν。
[0053] 步驟1 :設(shè)初始值C。: = 0, C中將要計(jì)算變量塊的初始位置s : = 0。對(duì)于目標(biāo)函 數(shù)（1)，設(shè) R : = F，res : = RTWR，RES := res，A: = 0〇
[0054] 步驟2 :若s+k>N，則s : = N轉(zhuǎn)到步驟3 ;否則，根據(jù)基本PSF計(jì)算出A的第 s+1，. . .，s+k列組成的子矩陣，并記為2，而方是Γ的第s+1，. . .，s+k組成的子矩陣， Zi :=(ATWA + BTaByl(A WR-BT〇TC0) , r, :=2RrWMs-X^WA^s-2C0^aBX s^XjBaBXx, Λ := i? -?，η := + , (^是 C 中的第 i 個(gè)元素，（cs+1，…，cs+k)T:= (c S+1，…，（^)4^:=5 + 1^/21 + 1，f1:= f Jr1, RES := RES-ri，轉(zhuǎn)到步驟2。
[0055] 步驟3 :若s-k < 0，s : = 0,轉(zhuǎn)到步驟4 ;否則，根據(jù)基本PSF計(jì)算出A的第 s-k+1，. . .，s列組成的子矩陣，并記為:i，而；§是Γ的第s-k+l，. . .，s組成的子矩陣， Xi := (A'WA +BTaBy1 (AWR-Bt〇TC0), r,-21^W^sWAXs-ICcrT raBXs -X/'B aBXs, -v4Xs，（^是 C 中的第 i 個(gè)元素 （c A!，· · ·，cs)T: = (c s k+1，…，cs)T+Xs， s:= 「A：/2~|- I，A: = f Jr1，RES : = RES-Ir1，轉(zhuǎn)到步驟 3。
[0056] 步驟4 :若f\< ε，停止；否則f 1: = 〇,轉(zhuǎn)到步驟2。
[0057] 最后，目標(biāo)函數(shù)也可以用并行計(jì)算求解。圖3是并行求解流程圖，具體如下：
[0058] 步驟0 :給定ε > 〇,子矩陣塊尺度k，和向量C的長(zhǎng)度N。0 < k << N。
[0059] 步驟 1 :設(shè)全局變量的初始值 C。:= 0，R : = F，res : = RTWR，rd: = 〇，RES : = res， f1:=〇。把變量C中按照變量編號(hào)順序分為「W/Αφ夬。除了最后一±夬，每塊有k個(gè)變量。 記這些塊為5,/ = 1,2,...,「7V/^，令$^丨是這些矩陣塊的集合。
[0060] 步驟2 :在集合W中選取變量塊，這些變量塊在矩陣A中對(duì)應(yīng)的子矩陣的非0行號(hào) 兩兩之間都不相同。記這些變量塊的集合為A。
[0061] 步驟3 :對(duì)于每個(gè)^ e Δ，用一個(gè)計(jì)算核計(jì)算。具體如下：記]和互分別是矩陣A和 Γ中與G相對(duì)應(yīng)列組成的子矩陣塊。令 X5 ：= (I1WA +B1 aBy1 (ArWR-Br〇TC0) , rs :=2^WAXs -Xj~AW~AXS-IC^YaBXs -Xj~B aBXs, R-R-AXs, €:=€ + 尤，匕：=f^+bRES:= RES-rs。
[0062] 步驟4 #:= W Δ，若p為空集，轉(zhuǎn)到步驟6 ;否則，轉(zhuǎn)到步驟2。
[0063] 步驟5 :重新把變量C進(jìn)行劃分成新的「Λ7/Γ1 + 1變量塊，使得每一個(gè)具有k個(gè) 元素塊與上次塊劃分的所有具有k個(gè)的元素塊至少有個(gè)變量不同。記這些塊為 ^^_ = 1,2,.··,「Λ7Α:]，令是這些矩陣塊的集合。轉(zhuǎn)到步驟3。
[0064] 步驟6 :若f\< ε，停止；否則f 1: = 〇,轉(zhuǎn)到步驟5。
[0065] 在本發(fā)明實(shí)施中，由于k值較小，所以在計(jì)算步驟中所有的矩陣塊運(yùn)算都可以進(jìn) 行稀疏存儲(chǔ)和計(jì)算，從而能加快計(jì)算速度。以上所述的具體實(shí)施例，對(duì)本發(fā)明的目的、技術(shù) 方案進(jìn)行了進(jìn)一步詳細(xì)說明，所應(yīng)理解的是，以上所述僅為本發(fā)明的具體實(shí)施例而已，并不 用于限制本發(fā)明，凡在本發(fā)明的精神和原則之內(nèi)，所做的任何修改、等同替換、改進(jìn)等，均應(yīng) 包含在本發(fā)明的保護(hù)范圍之內(nèi)。
【主權(quán)項(xiàng)】
1. 一種用反卷積方法從光纖光譜的二維CCD圖像中抽取一維光譜的方法。其特點(diǎn)是： 根據(jù)已知的點(diǎn)源擴(kuò)散函數(shù)（簡(jiǎn)記：PSF)，通過計(jì)算目標(biāo)函數(shù)的最小值，來反卷積整個(gè)光纖光 譜圖像，從而獲得一維光譜。該方法包括：光纖光譜中心跡線的獲得，每條光纖光譜所有波 長(zhǎng)處PSF的獲得，目標(biāo)函數(shù)的構(gòu)造，噪聲的抑制方法，目標(biāo)函數(shù)的串行和并行求解方法。2. 根據(jù)權(quán)利要求1所述的光纖光譜中心跡線的獲得方法，其特征在于，從較強(qiáng)的平場(chǎng) 定標(biāo)燈的二維光譜，利用重心法，求出每一條光纖光譜的每一行的重心位置，然后把每一行 重心用低階多項(xiàng)式擬合，得到每一條光纖光譜的中心軌跡。3. 根據(jù)權(quán)利要求1所述的PSF獲取方法，其特征在于，對(duì)于波長(zhǎng)定標(biāo)燈的每一條光纖光 譜圖像，首先獲取單根發(fā)射線的離散輪廓。對(duì)于發(fā)射線比較弱的情況，可以延長(zhǎng)曝光來提高 該發(fā)射線的信噪比。然后對(duì)離散的輪廓進(jìn)行歸一化，就獲得了該發(fā)射線在它所處CCD位置 處的離散的PSF。另外，也可以通過激光梳來獲得，只不過這種方法比較昂貴。如果想使用 光滑的PSF，可以通過利用B樣條曲面插值離散的PSF來獲得。我們把通過這種方法獲得的 PSF稱為基本的PSF。 對(duì)于每一條光纖光譜圖像，在沒有燈譜發(fā)射線的位置，我們可以通過與它相鄰的兩個(gè) 基本PSF線性插值獲得。 這些PSF都位于每一條光纖光譜的中心跡線上。4. 根據(jù)權(quán)利要求1所述的目標(biāo)函數(shù)的構(gòu)造方法，其特征在于，把圖像的噪聲分為高斯 噪聲和泊松噪聲兩種情況來構(gòu)造目標(biāo)函數(shù)。 高斯噪聲情況：其中，Nx表示圖像總列數(shù)；Ny表示圖像總行數(shù)；Nf表示圖像記錄的光纖光譜總數(shù)；F(k， m)表示圖像上位于第k行m列的CCD像素的計(jì)數(shù)；PSFuG^m)表示位于第i條光纖光譜在 第j行處的PSF在CCD上第k行第m列像素上的計(jì)數(shù)； Cl, j表示第i條光纖光譜在第j行 處的流量值，該值是要求解的值。 高斯噪聲下，目標(biāo)函數(shù)的矩陣形式寫法： - AC)T (F - AC) C 泊松噪聲下，目標(biāo)函數(shù)的矩陣形式寫法： mjn(F - ACf W(F -AC) C ' 其中，F(xiàn)是NxXNy行的列向量，F(xiàn)(k，m)是它的第kXNx+m個(gè)元素；A是N xXNy行，NfXNy 列的矩陣，？3匕#，111)是它的第1^凡+111行，第1\%+」列元素；(：是一個(gè)隊(duì)\%行的列向 量，h i是它第iXNy+j個(gè)元素；W是一個(gè)NxXNyR NxXNy列的對(duì)角矩陣，它的第kXNx+m行， kXNx+m列元素是F(k，m)，其它元素為0。5. 根據(jù)權(quán)利要求1所述的噪聲抑制方法，其特征在于，利用Tikhonov正則項(xiàng)來抑制噪 聲。可以使用所求變量的〇、1、2次等導(dǎo)數(shù)，及這些導(dǎo)數(shù)的組合來抑制噪聲。矩陣形式的目 標(biāo)函數(shù)如下： 高斯噪聲下，目標(biāo)函數(shù)為：其中，Γ是Tikhonov矩陣，源自導(dǎo)數(shù)的離散化；α是權(quán)重對(duì)角矩陣；Tikhonov項(xiàng)（rc) Ta (rc)對(duì)噪聲有抑制作用。 上面的公式的推導(dǎo)，都假設(shè)噪聲是獨(dú)立同分布的。若噪聲情況復(fù)雜，則刻畫噪聲的對(duì)稱 正定矩陣W會(huì)更復(fù)雜。無論何種情況，目標(biāo)函數(shù)可以寫為如下統(tǒng)一形式：(1) 其中，Γ是Tikhonov矩陣，源自導(dǎo)數(shù)的離散化；α是權(quán)重對(duì)角矩陣；Tikhonov項(xiàng)（rc) Ta (rc)對(duì)噪聲有抑制作用。6. 根據(jù)權(quán)利要求1所述的目標(biāo)函數(shù)的串行和并行求解方法，其特征在于，利用矩陣分 塊和迭代的方式求解目標(biāo)函數(shù)（1)，由于都是稀疏矩陣，矩陣在分塊計(jì)算過程中，利用了其 系數(shù)特性進(jìn)行存儲(chǔ)和計(jì)算。7. 根據(jù)權(quán)利要求6所述的目標(biāo)函數(shù)的串行求解方法，具體的計(jì)算步驟如下： 步驟〇:給定ε >〇,子矩陣塊尺度k，和向量C的長(zhǎng)度N。0<k<<N。 步驟1 :設(shè)初始值C。: = 0, C中將要計(jì)算的變量塊的初始位置S : = 0。對(duì)于目標(biāo)函數(shù) (1)，設(shè) R : = F，res : = RTWR，RES : = res，= 0〇 步驟2 :若s+k > N，則s : = N轉(zhuǎn)到步驟3 ;否則，記:i和云分別是矩陣A和Γ的第 s+1，…，s+k 列組成的子矩陣，尤:=(7#] + 礦6^)_1(7 呢-FarC0)，T?:=i? -IXi, 6 :=2Λ7詼IXi ，(^是 C 中的第 i 個(gè)元素， (cs+1，…，cs+k)T:= (c S+1，…，CsJkXvr=Hp/]"! + !，f1:= f Jr1, RES := RES-IT1，轉(zhuǎn)到 步驟2。 步驟3 :若s-k < 0，s : = 0,轉(zhuǎn)到步驟4 ;否則，記j和5分別是矩陣A和Γ的第 s-k+1，· · ·，s 列組成的子矩陣，Xs :=(7r^ +礦咖-7arC0)，i?:=i?-]iXs， η := - (2(^Τ?Χ5 +X/安《萬尤），' (^是 C 中的第 i 個(gè)元素， (cs k+1，…，cs)T: = (c s k+1，…，cs) T+Xs，S := S -「A:/2"| -1，A: = f Jr1，RES : = RES-IT1，轉(zhuǎn)到 步驟3。 步驟4 :若f\< ε，停止；否則f 1: = 〇,轉(zhuǎn)到步驟2。8.根據(jù)權(quán)利要求6所述的目標(biāo)函數(shù)的并行求解方法，其特征在于，在計(jì)算過程中，把要 求解的變量向量分成多個(gè)變量塊，每次并行計(jì)算那些沒有交叉污染的變量塊，直到算法收 斂到最小點(diǎn)。具體計(jì)算步驟如下： 步驟〇:給定ε >〇,子矩陣塊尺度k，和向量C的長(zhǎng)度N。0<k<<N。 步驟 1 :設(shè)全局變量的初始值 C。: = 0, R : = F，res : = RTWR，RES : = res，f1: = 0。把 變量C中按照變量編號(hào)順序分為塊。除了最后一塊，每塊有k個(gè)變量。記這些塊為 ，· = 1,2,…,「tv / q，令p :=后}是這些矩陣塊的集合。 步驟2 :在集合爐中選取變量塊，這些變量塊在矩陣A中對(duì)應(yīng)的子矩陣的非0行號(hào)兩兩 之間都不相同。記這些變量塊的集合為A。 步驟3 :對(duì)于每個(gè)Cs e Δ，用一個(gè)計(jì)算核計(jì)算。具體如下：記 j和云分別是矩陣A和Γ中與$相對(duì)應(yīng)列組成的子矩陣塊。令 Xi ：= (a W^A+ 1ToB)-1 (AWR-BtaTC0) ， R-=R-AXs ， Q :=Q+XS ， rs :=2RrWAXs-X^~AW~AXS-2C0 rTTaBXs-XsTBTaBX s, f1：= f !+r,, RES ：= RES-rs〇 步驟4少:=<?ΛΔ，若p為空集，轉(zhuǎn)到步驟6 ;否則，轉(zhuǎn)到步驟2。 步驟5 :重新把變量C進(jìn)行劃分成新的變量塊，使得每一個(gè)具有k個(gè)元素 塊與上次塊劃分的所有具有k個(gè)的元素塊至少有「&/2]-1個(gè)變量不同。記這些塊為 ^，/ = l,2,...,「iV/q，令是這些矩陣塊的集合。轉(zhuǎn)到步驟2。 步驟6 :若f\< ε，停止；否則f 1: = 〇,轉(zhuǎn)到步驟5。
【文檔編號(hào)】G01J3/28GK105890757SQ201410765824
【公開日】2016年8月24日
【申請(qǐng)日】2014年12月15日
【發(fā)明人】李廣偉, 張昊彤, 董義喬, 袁海龍, 雷亞娟, 白仲瑞, 楊卉沁
【申請(qǐng)人】中國(guó)科學(xué)院國(guó)家天文臺(tái)