基於貝葉斯線性估計的遙感圖像融合方法
2023-05-27 01:02:51 1
專利名稱：基於貝葉斯線性估計的遙感圖像融合方法
技術領域：
本發明屬於遙感圖像處理技術領域，具體涉及一種基於貝葉斯線性估計的遙感圖像融合方法。
背景技術：
由於遙感傳感器設計的限制，遙感圖像一般在空間解析度和光譜解析度之間進行折衷，具有高的光譜解析度的圖像一般不具備高的空間解析度，反之亦然。例如，Landsat ETM+傳感器提供的就是6幅30m空間解析度的多光譜波段圖像和一幅15m空間解析度的全色波段圖像。在實際應用中，那些同時具有高的空間解析度和高的光譜解析度的圖像能夠有效地提高解譯和分類的精度，因此不同解析度的遙感圖像的融合成為研究的熱點，尤其是低解析度的多光譜圖像和高解析度的全色波段圖像的融合。一般來說，融合後的圖像既要求保留多光譜圖像的光譜特性，又要求融入全色圖像的空間信息。
目前常見的融合方法有HIS方法[1]、PCA方法[2，3]以及小波變換方法[4，5]。HIS方法和PCA方法常常顯著地改變原多光譜圖像的光譜特性，而小波變換方法對於分解層次和小波基的選擇比較敏感，並且會因操作人員的不同，而有不同的融合效果。
近年來，基於統計參數估計的方法[6-9]開始應用於遙感圖像的融合中。Nishii等[8]假設高解析度多光譜圖像和全色圖像的概率分布服從聯合高斯分布，用條件均值作為估計。Hardie等[9]在上述聯合高斯分布假設的基礎上引入了高解析度多光譜圖像和低解析度多光譜圖像之間的觀測模型，得到了最大後驗概率(MAP)意義上的估計。Nishii方法和Hardie方法在低解析度多光譜圖像和全色圖像相關性不強的時候，難以融入全色圖像的空間信息。
針對以上問題，在遙感圖像融合的研究中，在增強空間細節的同時很好地保持光譜特性，以及保證算法的魯棒性成為目前研究的熱點。
發明內容
本發明的目的是提出一種基於貝葉斯線性估計的遙感圖像融合方法，以解決傳統的基於統計參數估計方法依賴於多光譜圖像和全色圖像相關係數的問題，增強空間細節，保持光譜特性。
本發明提出的基於貝葉斯線性估計的遙感圖像融合方法，具體步驟如下引入高解析度多光譜圖像和低解析度多光譜圖像之間的觀測模型，以及高解析度多光譜圖像和全色圖像之間的觀測模型，並將上述兩個觀測模型聯立成一個貝葉斯線性模型；應用貝葉斯高斯-馬爾科夫定理，計算得到線性最小均方誤差(LMMSE)意義上的高解析度多光譜圖像的估計。
下面對各步驟作進一步具體描述。
1、引入觀測模型假設同一地區分別被低解析度的多光譜波段傳感器和高解析度的全色波段傳感器拍攝到，所謂圖像解析度是指圖像中每個像素覆蓋地面的範圍，圖像解析度的高低是個相對的概念，根據需要可具體劃分。全色圖像按照如下方式排列成一維的列矢量x＝[x1，x2，…，x1，…，xN]T(1)其中xi表示全色圖像在空間位置i上的像素值，而N是全色圖像的像素數目。低解析度的多光譜圖像也按照類似的方式排列成一維的列矢量y=[y1T,y2T,,yjT,,yMT]T---(2)]]>其中yj表示低解析度的多光譜圖像在空間位置j上的像素值(有K個波段，yj＝[yj，1，yj，2，…，yj，k]T)，而M是低解析度多光譜圖像的像素數目。
假設高解析度的多光譜圖像存在，那麼它應該既包含多光譜圖像的光譜信息，又具有和全色圖像一樣的空間解析度，這裡用如下的一維列矢量來表示z=[z1T,z2T,,ziT,,zNT]T---(3)]]>其中zi表示高解析度的多光譜圖像在空間位置i上的像素值(有K個波段，zi＝[zi，1，zi，2，…，zi，k]T)，N是高解析度多光譜圖像的像素數目。
一般，低解析度的多光譜圖像可以認為是由高解析度的多光譜圖像(假如存在)通過低通濾波和降採樣過程得到的[10]。本發明引入高解析度的多光譜圖像與低解析度多光譜圖像之間的觀測模型，具體如下所示y＝Hz+u(4)其中u是隨機噪聲，它的均值是零，協方差矩陣是Cu，且與z是不相關的；H矩陣表示低通濾波和降採樣過程。
另外，我們在全色圖像和高解析度多光譜圖像之間，引入如下所示的觀測模型x＝GTz+v (5)其中v是隨機噪聲，它的均值為零，協方差矩陣是Cv，且與z是不相關的；G矩陣表示對高解析度多光譜圖像的K個波段做加權平均，權重因子如下所示
gl=ccl/l=1K|ccl|---(6)]]>其中ccl表示全色圖像與低解析度多光譜圖像的第l波段的相關係數。
由於式(4)和(5)提供的方程數目M×K+N小於待估計的未知量N×K(一般M＜N)，因此不能直接解出z。本發明將從貝葉斯線性模型的角度對上述兩個觀測模型進行推導，從而得到高解析度多光譜圖像z的估計量。
在估計高解析度的多光譜圖像時，觀測模型(4)和(5)可以按照如下的形式聯立成一個貝葉斯線性模型yx=HGTz+uv---(11)]]>2、應用貝葉斯高斯-馬爾科夫定理，進行線性估計假設數據由貝葉斯線性模型來描述＝Aθ+w(7)其中是L×1數據矢量，A是已知的L×p觀測矩陣，θ是p×1的隨機矢量。θ的現實是要估計的，它的均值是E(θ)，協方差矩陣是Cθ。w是L×1的隨機矢量，它的均值是零，協方差矩陣是Cw，且與θ是不相關的。
首先假定θ的估計量 可通過數據集按下式求得， 選擇加權係數aj和bi來使貝葉斯均方誤差(Bayesian mean square error，BMSE)BMSE(^)=E[(-^)2]---(9)]]>最小，得到的估計量稱為線性最小均方誤差(Linear minimum mean square error，LMMSE)估計量(貝葉斯高斯-馬爾可夫定理[11])，如下所示 此時BMSE(^l)=[(C0+ATCw-1A)-1]ii.]]>由於通常情況下θ不能完美地表示為j的線性組合，所以LMMSE估計量不是最佳的，但是在實際中它是相當有用的，因為它具有閉合形式的解，並且只與均值和協方差有關。
對於貝葉斯線性模型(11)，應用(10)式中的LMMSE估計量可以得到高解析度的多光譜圖像的估計(這裡，令n＝[uT，vT]T)
z^=E(z)+CzHGTT(HGTGzHGTT+Gn)-1y-E(y)x-E(x)---(12)]]>其中的統計參數估計如下在式(12)中要得到估計量 需要知道高解析度多光譜圖像的均值E(z)和協方差Cz。這裡假設高解析度多光譜圖像的像素之間是相互獨立的，因此整個圖像的均值和協方差可以由各個像素的均值和協方差組成[12]，具體如下所示E(z)=[E(z1)T,E(z2)T,,E(zi)T,,E(zN)T]---(13)]]> 其中E(z)是使用雙線性插值(B)後的多光譜圖像，具體為E(z)＝B(y) (15)為了估計Cz，我們用矢量量化算法對上述矢量E(z1)根據歐式距離進行分類，計算每一類矢量集的協方差矩陣，並把它作為類中各矢量對應的協方差矩陣。
另外，式(12)中的E(y)的估計是通過對插值圖像低通濾波和降採樣(H)得到的，E(y)＝H(E(z))＝H(B(y)) (16)而E(x)的估計是通過對全色波段圖像進行低通濾波和降採樣，然後再雙線性插值得到的，E(x)＝B(H(x))(17)在實際的計算中，如果協方差矩陣維數比較大，(12)式中的求逆運算將會產生困難，因此，在這種情況下，本發明將圖像分割成許多圖像小塊再進行估計。
上述的貝葉斯線性估計方法對於融合前圖像的波段數目沒有限制，低解析度的多光譜圖像y的波段數目可以大於3，而高解析度的全色圖像既可以是單波段也可以是多波段的。因此，貝葉斯線性估計方法既適用於多波段的多光譜圖像和單波段的全色圖像的融合，也適用於多波段的超光譜圖像和多波段的多光譜圖像的融合。
圖1為基於貝葉斯線性估計的遙感圖像融合方法的幾何解釋。
圖2為Landsat ETM+的全色圖像和多光譜圖像融合實驗結果。其中，圖2(a)是30m空間解析度的多光譜圖像，圖2(b)是30m空間解析度的全色圖像，圖2(c)是120m空間解析度的多光譜圖像，圖2(d)是HIS方法，圖2(e)是PCA方法，圖2(f)是小波變換方法，圖2(g)是Nishii方法，圖2(h)是Hardie方法，圖2(i)是貝葉斯方法圖3為SPOT的全色圖像和TM的多光譜圖像融合實驗結果。其中，圖3(a)是30m空間解析度的多光譜圖像，圖3(b)是30m空間解析度的全色圖像，圖3(c)是150m空間解析度的多光譜圖像，圖3(d)是HIS方法，圖3(e)是PCA方法，圖3(f)是小波變換方法，圖3(g)是Nishii方法，圖3(h)是Hardie方法，圖3(i)是貝葉斯方法具體實施方式
以下通過實施例，進一步對發明中的各個組成加以描述。
1設置觀測模型低解析度的多光譜圖像可以認為是由高解析度的多光譜圖像(假如存在)通過低通濾波和降採樣過程得到的[10]，該觀測模型如下所示y＝Hz+u(4)其中u是隨機噪聲，它的均值是零，協方差矩陣是Cu，且與z是不相關的；H矩陣表示低通濾波和降採樣過程。
在全色圖像和高解析度多光譜圖像之間，存在著如下所示的觀測模型x＝GTz+v (5)其中v是隨機噪聲，它的均值為零，協方差矩陣是Cv，且與z是不相關的；G矩陣表示對高解析度多光譜圖像的K個波段做加權平均，權重因子如下所示gl=ccl/l=1K|ccl|---(6)]]>其中ccl表示全色圖像與低解析度多光譜圖像的第l波段的相關係數。
2將觀測模型聯繫成貝葉斯線性模型在估計高解析度的多光譜圖像時，觀測模型(4)和(5)可以按照如下的形式聯立成一個貝葉斯線性模型yx=HGTz+uv---(11)]]>3應用貝葉斯線性估計對於一個貝葉斯線性模型(11)，應用(10)式中的LMMSE估計量可以得到高解析度的多光譜圖像的估計(這裡，令n＝[uT，vT]T)z^=E(z)+CzHGTT(HGTGzHGTT+Gn)-1y-E(y)x-E(x)---(12)]]>4估計統計參數假設高解析度多光譜圖像的像素之間是相互獨立的，因此整個圖像的均值和協方差可以由各個像素的均值和協方差組成[12]，如下所示E(z)＝[E(z1)T，E(z2)T，…，E(zi)T，…，E(zN)T] (13) 其中E(z)是使用雙線性插值(B)後的多光譜圖像，E(z)＝B(y)(15)為了估計Cz，我們用矢量量化算法對上述矢量E(zi)根據歐式距離進行分類，計算每一類矢量集的協方差矩陣，並把它作為類中各矢量對應的協方差矩陣。
另外，式(12)中的E(y)的估計是通過對插值圖像低通濾波和降採樣(H)得到的，E(y)＝H(E(z))＝H(B(y)) (16)而E(x)的估計是通過對全色波段圖像進行低通濾波和降採樣，然後再雙線性插值得到的，如下所示E(x)＝B(H(x)) (17)對本發明的上述方法，進行了仿真計算。具體的仿真條件如下(1)Landsat 7 ETM+的全色圖像和多光譜圖像；(2)SPOT的全色圖像和TM的多光譜圖像。
圖2顯示的是Landsat 7 ETM+傳感器在2000年6月14日拍攝的上海地區的多光譜圖像和全色圖像。其中，全色圖像具有15m的空間解析度，而多光譜圖像具有30m的空間解析度，這裡以第3，2和1波段分別作為RGB通道。
由於Landsat 7 ETM+不提供15m解析度的真實多光譜圖像用來比較，因此我們很難評價各種方法得到的融合效果。為了解決這個問題，我們將全色圖像和多光譜圖像的空間解析度分別退化到30m和120m。通過對上述退化圖像進行融合，並將融合結果與原30m解析度的多光譜圖像進行比較。
圖3顯示了在1995年10月26日拍攝的Hanoi地區的SPOT衛星的全色圖像和TM的多光譜圖像(http//earth.esa.int/mcities/images/cases)。本發明將SPOT的全色圖像從10m的解析度退化到30m，TM的多光譜圖像從30m的解析度退化到150m。這裡以TM的第3，2，和1波段分別作為RGB通道。
本發明所提出的方法將和以下方法做比較(1)傳統的HIS方法、PCA方法；
(2)小波變換方法，選用4階Daubechies小波基，分解的層次是3層。
(3)Nishii方法和Hardie方法。
實驗結果如下1、Landsat 7 ETM+的全色圖像和多光譜圖像從視覺效果上分析，圖2中基於HIS方法和PCA方法的融合圖像顯著地改變了真實多光譜圖像的光譜特性，其融合效果明顯劣於其它圖像融合方法，因此在以下的基於統計參數的定量比較中，不再考慮HIS方法和PCA方法。
小波變換方法導致右邊的水體中出現了一些斑點。Nishii方法和Hardie方法的融合圖像很模糊。這裡需要指出的是低解析度多光譜圖像和全色圖像的相關係數分別是0.57、0.18和-0.20。貝葉斯方法的融合圖像比較接近於真實的多光譜圖像，不僅增強了空間細節，而且很好地保持了原多光譜圖像的光譜特性。
以下部分我們定量地比較小波變換方法、Nishii方法和Hardie方法以及貝葉斯方法在增強空間細節和保持光譜特性上的統計參數變化情況。
為了衡量上述方法增強空間細節的能力，我們計算融合圖像減去120m解析度的多光譜圖像得到的差值圖像在各個波段上的標準差(SDD)，如表1所示。其中第一列顯示了真實的多光譜圖像(30m解析度)的SDD參數。在表1中，Nishii方法和Hardie方法在各通道上的SDD參數都遠遠小於真實圖像的值，說明Nishii方法和Hardie方法在本實驗中不能有效地增強空間細節。小波變換方法在各通道上的SDD參數是相同的，這說明小波變換方法主要取決於分解層次，而對各通道上的具體情況不作針對性的處理。貝葉斯方法的SDD參數在各通道上比較接近於真實圖像的值，這說明通過高解析度多光譜圖像的協方差估計來決定融合圖像在各通道上的幅度的做法是合理的。
表1 各種方法增強空間細節的統計參數(Landsat ETM+)
為了衡量上述方法保持光譜特性的能力，我們採用如下的統計參數(1)峰值信噪比[9](PSNR)用來衡量圖像灰度的峰值和兩幅圖像的誤差之間的比率，定義如下PSNR＝20×log10(b/rms)(20)其中b表示圖像灰度的峰值，rms是兩幅圖像的誤差均方根，峰值信噪比的單位是db。一般峰值信噪比越大，兩幅圖像之間的差異越小。
(2)相關係數(CC)的定義如下C(f,g)=[(f(i,j)-f)(g(i,j)-g)][(f(i,j)-f)2][(g(i,j-g))2]---(21)]]>其中f(i，j)和g(i，j)表示圖像的灰度，f和g是圖像的均值。一般相關係數越高，兩幅圖像越相似。峰值信噪比和相關係數在融合圖像和真實的多光譜圖像的各波段上分別計算。
(3)相對全局誤差[5](ERGAS)如下所示，RASE參數越低，光譜的扭曲程度越小。
ERGAS=100hl1Kk=1Krms(K)/mz(K)---(22)]]>其中，l和h是融合前和融合後多光譜圖像的解析度(這裡分別取120m和30m)，rms(K)表示融合圖像和真實的多光譜圖像在各波段上的誤差均方根，mz(K)是真實的多光譜圖像在各波段上的均值。一般相對全局誤差越小，說明融合圖像的光譜特性保持的越好。
上述統計參數如表2所示。貝葉斯方法在各通道上都具有最高的峰值信噪比，和最大的相關係數，而相對全局誤差是最小的，說明貝葉斯方法的融合圖像和真實的多光譜圖像之間的差異最小，因此貝葉斯方法在光譜特性的保持上優於小波變換方法、Nishii方法和Hardie方法。
表2 各種方法保持光譜特性的統計參數(Landsat ETM+)
2、SPOT的全色圖像和TM的多光譜圖像從視覺效果上分析，圖3中的HIS方法和PCA方法顯著地改變了多光譜圖像的光譜特性，尤其是PCA方法將河流的色彩轉換。Nishii方法和Hardie方法的融合圖像比較模糊。這裡，低解析度多光譜圖像和全色圖像的相關係數分別是0.52、0.35和0.29。小波變換方法和貝葉斯方法的融合效果相對比較好。
在表3中，Nishii方法和Hardie方法的SDD參數都很小，說明它們對於空間解析度的提高確實很有限。小波變換方法在各通道上具有相同的SDD參數，仍然沒有對不同通道上的具體情況作有針對性的處理。貝葉斯方法在各通道上的SDD參數和真實多光譜圖像的值幾乎相同。
另外，在表4中，貝葉斯方法在各通道上都具有最高的峰值信噪比和最大的相關係數，而相對全局誤差是最小的，說明貝葉斯方法的融合圖像和真實的多光譜圖像之間的差異較小，因此，可以說貝葉斯方法的性能優於小波變換方法、Nishii方法和Hardie方法。
表3 各種方法增強空間細節的統計參數(SPOT和TM)
表4 各種方法保持光譜特性的統計參數(SPOT和TM)
綜上所述，貝葉斯方法解決了Nishii方法和Hardie方法依賴於多光譜圖像和全色圖像的相關係數的問題，同時實驗結果證明了貝葉斯方法的融合性能明顯優於HIS方法、PCA方法和小波變換方法。
最後，特別值得指出的是，一些方法如小波變換方法，由於參數設置和操作人員的不同，會導致融合結果不同，而本發明所提出的方法將自動設置參數，在不需要人為幹預的情況下仍能取得較好的融合效果。
參考文獻[1]J W Carper，T M Lillesand，R W Kiefer.The use of intensity-hue-saturation transformationfor merging SPOT panchromatic and multispectral image data[J].PhotogrammetricEngineering and Remote Sensing，1990，56459-467. V K Shettigara.A generalized component substitution technique for spatial enhancement ofmultispectral images using a higher resolution data set[J].Photogrammetric Engineering andRemote Sensing，1992，58561-567. P S Chavez，S C Sides，J A Anderson.Comparison of three different methods to mergemulti-resolution and multispectral daraLandsat TM and SPOT panchromatic[J].Photogrammetric Engineering and Remote Sensing，1991，57(3)295-303. J X Otazu，O Fors，et al.Multiresolution-based image fusion with additive waveletdecomposition[J].IEEE Transactions on Geoscience and Remote Sensing，1999，371204-1211. M A González，J L Saleta，R G Catalán，et al.Fusion of mnltispectral and panchromaticimages using improved IHS and PCA mergers based on wavelet decomposition[J].IEEETransactions on Geoscience and Remote Sensing，2004，23(18)1291-1299. J C Price.Combining panchromatic and multispectral imagery from dual resolution satelliteinstruments[J].Remote Sensing of Environment，1987，21119-128. C K Munechika，J S Warnick，C Salvaggio，et al.Resolution enhancement of multispectralimage data to improve classification accuracy[J].Photogrammetric Engineering and RemoteSensing，1993，5967-72. R Nishii，S Kusanobu，S Tanaka.Enhancement of low resolution image based on highresolution bands[J].IEEE Transactions on Geoscience and Remote Sensing，1996，341151-1158. R C Hardie，M T Eismann，G L Wilson.MAP estimation for hyperspectral image resolutionenhancement using an auxiliary sensor.IEEE Transactions on image processing，2004，13(9)1174-1184. M T Eismann，R C Hardie.Application of stochastic mixing model to hyperspectralresolution enhancement[J].IEEE Transactions on Geoscience and Remote Sensing，2004，42(9)1924-1933. S M Kay，Fundamentals of statistical signal processingEstimation theory[M].EnglewoodsCliffs，NJPrentice-Hall，1993391-392. H Eves.Elementary matrix theory[M].New YorkDover，1966107. J Zhou，D L Civco，J A Silander.A wavelet transform method to merge Landsat TM andSPOT panchromatic data[J].International Journal of Remote Sens.，1998，19(4)743-757.
權利要求
1.一種基於貝葉斯線性估計的遙感圖像融合方法，其特徵在於引入高解析度多光譜圖像和低解析度多光譜圖像之間的觀測模型，以及高解析度多光譜圖像和全色圖像之間的觀測模型，並將上述兩個觀測模型聯立成一個貝葉斯線性模型；應用貝葉斯高斯-馬爾科夫定理，計算得到線性最小均方誤差LMMSE意義上的高解析度多光譜圖像的估計。
2.根據權利要求1所述的基於貝葉斯線性估計的遙感圖像融合方法，其特徵在於所述的高解析度多光譜圖像和低解析度多光譜圖像之間的觀測模型如下y＝Hz+u (4)其中u是隨機噪聲，它的均值是零，協方差矩陣是Cu，且與z是不相關的；H矩陣表示低通濾波和降採樣過程；y為低解析度的多光譜圖像排成一維的列矢量y=[y1T,y2T,...,yjT,...,yMT]T---(2)]]>其中yj表示低解析度的多光譜圖像在空間位置j上的像素值(有K個波段，yj＝[yj，1，yj，2，…，yj，K]T)，而M是低解析度多光譜圖像的像素數目；z為高解析度的多光譜圖像排成一維的列矢量z=[z1T,z2T,...,ziT,...,zNT]T---(3)]]>其中zi表示高解析度的多光譜圖像在空間位置i上的像素值(有K個波段，zi＝[zi，1，zi，2，…，zi，K]T)，N是高解析度多光譜圖像的像素數目；所述高解析度多光譜圖像和全色圖像之間的觀測模型如下x＝GTz+v(5)其中v是隨機噪聲，它的均值為零，協方差矩陣是Cv，且與z是不相關的；G矩陣表示對高解析度多光譜圖像的K個波段做加權平均，權重因子如下所示gl=ccl/l=1K|ccl|---(6)]]>其中ccl表示全色圖像與低解析度多光譜圖像的第l波段的相關係數；x為全色圖像排成的一維矢量x＝[x1，x2，…，xi，…，xN]T(1)其中xi表示全色圖像在空間位置i上的像素值，而N是全色圖像的像素數目；觀測模型(4)和(5)按照如下的形式聯立成一個貝葉斯線性模型yx=HGTz+uv---(11)]]>
3.根據權利要求2所述的基於貝葉斯線性估計的遙感圖像融合方法，其特徵在於高解析度多光譜圖像的估計式如下z^=E(z)+CzHGTTHGTCzHGTT+Cn-1y-E(y)x-E(x)---(12)]]>其中，E(z)和Cz分別為高解析度多光譜圖像的均值的協方差，它們由各個像素的均值和協方差組成，具體如下所示E(z)＝[E(z1)T，E(z2)T，…，E(zi)T，…，E(zN)T](13) 其中E(z)是使用雙線性插值(B)後的多光譜圖像，具體為E(z)＝B(y)(15)為了估計Cz，用矢量量化算法對上述矢量E(zi)根據歐式距離進行分類，計算每一類矢量集的協方差矩陣，並把它作為類中各矢量對應的協方差矩陣；式(12)中的E(y)的估計是通過對插值圖像低通濾波和降採樣(H)得到的E(y)＝H(E(z))＝H(B(y)) (16)而E(x)的估計是通過對全色波段圖像進行低通濾波和降採樣，然後再雙線性插值得到的E(x)＝B(H(x)) (17)。
全文摘要
本發明屬於遙感圖像處理技術領域，具體為一種基於統計參數估計的遙感圖像融合方法。該方法引入高解析度多光譜圖像和低解析度多光譜圖像之間的觀測模型，以及高解析度多光譜圖像和全色圖像之間的觀測模型，並將上述兩個觀測模型聯立成一個貝葉斯線性模型。通過應用貝葉斯高斯—馬爾科夫定理，得到線性最小均方誤差意義上的高解析度多光譜圖像的估計。本發明不僅能夠增強空間細節，而且很好地保持了光譜特性，其性能優於傳統的HIS方法、PCA方法和小波變換方法以及現存的基於統計參數估計的Nishii方法和Hardie方法。新方法可為改善遙感圖像的目視判讀精度，提高信息清晰度和可靠性上提供新的有效的技術支持。
文檔編號G01S7/48GK1808181SQ20061002411
公開日2006年7月26日 申請日期2006年2月23日 優先權日2006年2月23日
發明者葛志榮, 王斌, 張立明 申請人:復旦大學