圖像特征提取方法

1、. . . . 圖像特征提取方法摘 要  特征提取是計算機視覺和圖像處理中的一個概念。它指的是使用計算機提取圖像信息，決定每個圖像的點是否屬于一個圖像特征。特征提取的結果是把圖像上的點分為不同的子集，這些子集往往屬于孤立的點、連續(xù)的曲線或者連續(xù)的區(qū)域。 至今為止特征沒有萬能和精確的圖像特征定義。特征的精確定義往往由問題或者應用類型決定。特征是一個數(shù)字圖像中“有趣”的部分，它是許多計算機圖像分析算法的起點。因此一個算法是否成功往往由它使用和定義的特征決定。因此特征提取最重要的一個特性是“可重復性”：同一場景的不同圖像所提取的特征應該是一樣的。    

2、60;    特征提取是圖象處理中的一個初級運算，也就是說它是對一個圖像進行的第一個運算處理。它檢查每個像素來確定該像素是否代表一個特征。假如它是一個更大的算法的一部分，那么這個算法一般只檢查圖像的特征區(qū)域。作為特征提取的一個前提運算，輸入圖像一般通過高斯模糊核在尺度空間中被平滑。此后通過局部導數(shù)運算來計算圖像的一個或多個特征。 常用的圖像特征有顏色特征、紋理特征、形狀特征、空間關系特征。當光差圖像時，常?？吹降氖沁B續(xù)的紋理與灰度級相似的區(qū)域，他們相結合形成物體。但如果物體的尺寸很小或者對比度不高，通常要采用較高的分辨率觀察：如果物體的尺寸很大或對比度很強，只需要

3、降低分辨率。如果物體尺寸有大有小，或對比有強有弱的情況下同事存在，這時提取圖像的特征對進行圖像研究有優(yōu)勢。常用的特征提取方法有：Fourier變換法、窗口Fourier變換（Gabor)、小波變換法、最小二乘法、邊界方向直方圖法、基于Tamura紋理特征的紋理特征提取等。設計容課程設計的容與要求（包括原始數(shù)據(jù)、技術參數(shù)、條件、設計要求等）：一、課程設計的容本設計采用邊界方向直方圖法、基于PCA的圖像數(shù)據(jù)特征提取、基于Tamura紋理特征的紋理特征提取、顏色直方圖提取顏色特征等等四種方法設計。（1）邊界方向直方圖法由于單一特征不足以準確地描述圖像特征,提出了一種結合顏色特征和邊界方向特征的圖像檢

4、索方法.針對傳統(tǒng)顏色直方圖中圖像對所有像素具有一樣重要性的問題進行了改進,提出了像素加權的改進顏色直方圖方法;然后采用非分割圖像的邊界方向直方圖方法提取圖像的形狀特征,該方法相對分割方法具有簡單、有效等特點,并對圖像的縮放、旋轉以與視角具有不變性.為進一步提高圖像檢索的質(zhì)量引入相關反饋機制,動態(tài)調(diào)整兩幅圖像相似度中顏色特征和方向特征的權值系數(shù),并給出了相應的權值調(diào)整算法.實驗結果表明,上述方法明顯地優(yōu)于其它方法.小波理論和幾個其他課題相關。所有小波變換可以視為時域頻域的形式，所以和調(diào)和分析相關。所有實際有用的離散小波變換使用包含有限脈沖響應濾波器的濾波器段(filterbank)。構成CWT的

5、小波受海森堡的測不準原理制約，或者說，離散小波基可以在測不準原理的其他形式的上下文中考慮。通過邊緣檢測，把圖像分為邊緣區(qū)域和非邊緣區(qū)域，然后在邊緣區(qū)域進行邊緣定位根據(jù)局部區(qū)域邊緣的直線特性，求得小鄰域直線段的高精度位置；再根據(jù)邊緣區(qū)域邊緣的全局直線特性，用線段的中點來擬合整個直線邊緣，得到亞像素精度的圖像邊緣在擬合的過程中，根據(jù)直線段轉角的變化剔除了噪聲點，提高了定位精度并且，根據(jù)角度和距離區(qū)分出不同直線和它們的交點，給出了圖像精確的矢量化結果圖像的邊界是指其周圍像素灰度有階躍變化或屋頂變化的那些像素的集合，邊界廣泛的存在于物體和背景之間、物體和物體之間，它是圖像分割所依賴的重要特征邊界方向直

6、方圖具有尺度不變性，能夠比較好的描述圖像的大體形狀邊界直方圖一般是通過邊界算子提取邊界，得到邊界信息后，需要表征這些圖像的邊界，對于每一個邊界點，根據(jù)圖像中該點的梯度方向計算出該邊界點處法向量的方向角，將空間量化為M級，計算每個邊界點處法向量的方向角落在M級中的頻率，這樣便得到了邊界方向直方圖圖像中像素的梯度向量可以表示為 ( ，)，)， ( ，)，) ，其中Gx( ，)，)，G ( ，)，)可以用下面的Sobel算子8：( ，Y)=，( 一1，Y+1)+2×，( ，Y+1)+，( +1，Y+1)一，( 一1，Y一1)一2× ，( ，Y一1)一，( +1，Y一1)( ，Y)

7、=，( +1，Y一1)+2×，( +1，Y)+，( +1，Y+1)一，( 一1，Y一1)一2× ，( 一1，Y)一，( 一1，Y+1) (4)I(x，y)表示在圖像的(x，y)點像素的亮度。為了減少由于數(shù)字化過程中產(chǎn)生的噪聲的影響，邊界方向直方圖是基于局部像素梯度向量的平均值因為相反方向的梯度可能會相互抵消，所以并不能直接對局部鄰居像素點的梯度向量平均解決這個問題的一個方法就是在計算平均值前，先對用復數(shù)表示的向量進行平方運算，等價于把梯度向量的角度增加一倍角度增加一倍的相反方向的兩個梯度向量分別指向它們的對立梯度向量，從而互相得到增強；而處于正交的兩個梯度將會互相抵消，經(jīng)過

8、平均運算后變回到它們的單一角度的表示。平方的向量可以表示為：(G +_，G ) =G 一G +2_，G G (5)梯度平方的平均值可以通過對局部鄰居(利用非均勻窗口w)的平均計算得到，DFx：(G +G；)，DFy：(2GxGr) (6)W W現(xiàn)在，平均梯度方向 0，180可以由下面的公式得到：r l80， DFx=0 and DFx=0o(x，Y)=l80 DF (7) 麗， mM因此，邊界方向直方圖可以定義為：H(k)：H(k)+1，f ( ，Y)：k and k 0，180 (8)= 180是一個特例它并不意味著像素的邊界方向是水平方向且與 =0的一樣它表示在這個像素的周。圍沒有亮度變化

9、邊界方向直方圖之間距離D 的計算采用與上面顏色直方圖距離之間的一樣的計算公式(3)而得到。(2)基于PCA的圖像數(shù)據(jù)特征提取主成分分析 （ Principal Component Analysis ， PCA ）是一種掌握事物主要矛盾的統(tǒng)計分析方法，它可以從多元事物中解析出主要影響因素，揭示事物的本質(zhì)，簡化復雜的問題。計算主成分的目的是將高維數(shù)據(jù)投影到較低維空間。給定 n 個變量的 m 個觀察值，形成一個 n m 的數(shù)據(jù)矩陣， n 通常比較大。對于一個由多個變量描述的復雜事物，人們難以認識，那么是否可以抓住事物主要方面進行重點分析呢？如果事物的主要方面剛好體現(xiàn)在幾個主要變量上，我們只需要將這幾

10、個變量分離出來，進行詳細分析。但是，在一般情況下，并不能直接找出這樣的關鍵變量。這時我們可以用原有變量的線性組合來表示事物的主要方面， PCA 就是這樣一種分析方法。圖像的有效檢索手段逐漸引起人們的重視,傳統(tǒng)的采用基于關鍵字或描述性文本的數(shù)據(jù)庫檢索方式,已遠遠不能滿足人們的需要,近年來,國際上廣泛開展了基于容的圖像檢索(CBIR) 的研究,CBIR 是指根據(jù)圖像的顏色、形狀、紋理等特征以與這些特征的組合來查詢圖像,是計算機圖像處理和數(shù)據(jù)庫技術的有效結合。顏色是彩色圖像的最顯著特征之一。因此,基于顏色的查詢是基于容的圖像檢索中最基本的方法。這種查詢可針對任何類型的彩色圖像。目前,人們已提出了許多

11、種借助于顏色特征對圖像進行檢索的方法。這些方法常用的彩色空間有RGB 和HSV ,提取的主要特征是顏色直方圖。常用的計算直方圖之間距離的方法有直方圖相交法、距離法等，以直方圖作為圖像的顏色特征進行圖像檢索,需存儲大量的冗余特征信息,不能精確而簡潔地描述圖像。本文提出了一種基于HSV 顏色模型,利用主元分析法(PCA) 對圖像的顏色進行特征提取的方法,大大降低了特征向量的維數(shù),減小了特征數(shù)據(jù)庫的規(guī)模。圖像的顏色有多種表示方式,其中HSV 顏色模型是一種適合肉眼分辨的模型,它把彩色信息表示為三種屬性: 色調(diào)h 、飽和度s 和亮度v 。HSV 模型的色調(diào)h 是由顏色名稱來辨別的,如紅、橙、綠,它用角

12、度0°360°度量;亮度v 是顏色的明暗程度, 通常用百分比度量, 從黑0 到白100 %;飽和度s 指顏色的深淺,用百分比來度量,為從0 到完全飽和的100 %·這種顏色模型用Munsell 三維空間坐標系統(tǒng)表示,因坐標之間的心理感知獨立性,因此,可以獨立感知各顏色分量的變化;且這種顏色模型具有線性伸縮性,可感知的顏色差是與顏色分量的相應樣值上的歐氏距離成比例的·在CBIR 中應用這種模型更適合用戶的主元分析法(PCA) 對L 的降維處理由節(jié)1. 2 可得到L 矩陣, L 為一m ×n 矩陣, m 、n 由圖像尺寸決定。主元分析法( PCA)

13、 是一種線性降維技術,其基本思想是通過對數(shù)據(jù)協(xié)方差矩陣的分解,在其n 個特征值中取前a 個特征值, ( n - a) 個特征值被濾出。a 遠小于n ,而a個特征值對應的特征向量構成負荷矩陣P Rn ×a , L 到低維空間的投影就包含在得分矩陣中T = L P (3)由T 返回到n 維空間的投影為L = TPT (4)由L 成的子空間稱為得分空間,其所包含的信息量近似原空間L 中所包含的信息量9 。由矩陣L 可求出其協(xié)方差矩陣S 的特征值分解S = 1/ ( n - 1) L T L = VV T (5)包含幅值遞減的非負實特征值(1 2 n 0) ·為了最優(yōu)地獲取數(shù)據(jù)的變

14、化量,這里a =6 ,即取前6 個特征值。利用主元分析法(PCA) 對彩色圖像顏色特征進行提取,可將原圖像從m ×n ×3 降低至m ×a ( a 遠小于n) ,大大降低了特征向量的維數(shù),減小了特征數(shù)據(jù)庫的規(guī)模,為圖像的檢索建立了良好的基礎·將提取的顏色特征與紋理特征相結合作為支持向量機(SVM) 的輸入,同時引入相關反饋來對圖像進行分類與檢索,以便近一步提高檢索效率是有待進一步研究的容。(3)基于Tamura紋理特征的紋理特征提取隨著多媒體信息的廣泛應用,并在數(shù)據(jù)庫系統(tǒng)和計算機視覺兩大研究領域的共同推動下,圖像檢索技術己逐漸成為一個非?；钴S的研究領域。

15、如何有效地對這些圖像進行分析、存儲和檢索是一個亟待解決的問題?；谌莸膱D像檢索技術能有效的解決這一問題,成為研究的熱點?；谌莸膱D像檢索(CBIR)不同于傳統(tǒng)的基于文本的圖像檢索,它實際上是一種模糊查詢技術,通過對圖像提取一定的特征,找出在特征空間中與查詢要求接近的圖像,從而實現(xiàn)在圖像數(shù)據(jù)庫中自動地、智能地檢索、查詢和管理圖像。在CBIR系統(tǒng)中,特征的提取和匹配算法是決定圖像檢索結果的關鍵。在已經(jīng)存在的顏色特征,紋理特征,形狀特征等幾種特征提取方法中,由于紋理特征能夠描述圖像的平滑,稀疏,規(guī)則性等特性,因此本文決定采用紋理特征作為圖像檢索的依據(jù)。 首先,介紹了國外基于容的圖像檢索系統(tǒng)的研究現(xiàn)狀

16、,從基于容的圖像檢索系統(tǒng)的基本原理入手,論述了基于容的圖像檢索系統(tǒng)的通用框架、關鍵技術、檢索特點以與應用領域。 其次,描述了紋理和紋理特征,介紹了幾種常用的紋理分析方法,有灰度直方圖、邊緣方向直方圖、空間灰度共生矩陣、Tamura紋理特征、Gabor濾波器和Gabor小波以與圖像紋理譜分析方法,本文在傳統(tǒng)的紋理特征提取方法的基礎上,利用Gabor小波變換在數(shù)字圖像中局部區(qū)域的頻率,方向信息的優(yōu)異性能,和共生矩陣對圖像的整體區(qū)域有著較好的處理效果的特點,嘗試了結合灰度共生矩陣和Gabor濾波器來提取紋理特征的方法,并對提取的特征進行高斯歸一化處理。 另外,在圖像相似性度量方面,闡述了幾何矩陣模型

17、和集合理論模型,給出了兩種多特征組合的相似性度量結構,以與相似性度量的方法。同時介紹了在基于容的圖像檢索技術中相關性反饋技術。 最后,本文設計了一個基于紋理特征的圖像檢索原型系統(tǒng),并對系統(tǒng)的各個模塊與其功能進行了介紹,對本文提出的檢索技術進行了仿真實驗,給出了實驗結果,并用排序平價法和平均檢索率對算法進行了評價。Tamura紋理特征表示方法,并對其相似性度量方法進行了改進,使用直方圖表示了最重要的粗糙度、對比度、方向度特征。本文還將粗糙度直方圖改進為能夠適應圖片尺寸變化的方式,另外使用對數(shù)非均勻量化方式將對比度特征形成直方圖,該處理方法能夠更適應人眼的認知特點。同時,對每種紋理特征都給出了相似

18、性度量仿真實例,證實該方法具有較好效果。最后使用三種特征的加權距離對圖像的紋理特征相似性進行度量,并給出了圖像檢索的結果和性能評價。 3)興趣點特征提取:本文綜合了Harris興趣點標定方法和區(qū)域矩不變特性對圖像興趣點特征進行表達。該特征能夠適應光照、視角、縮放等場景的變化,文中給出圖片實例對這一特點進行了驗證。本文根據(jù)實驗提出了圖像的興趣點相似性度量方法,并給出了圖像檢索的結果和性能評價。 本文第二部分研究了區(qū)域相似性評價和圖片共性特征歸納方法。區(qū)域分割是該過程的基礎,本文介紹了一種簡單的基于顏色和紋理的區(qū)域分割法,并提出了圖像空間結構相似性度量方法。在對具有相似容的圖片進行共性特征提取過程

19、中,根據(jù)不同的應用場景選取了不同的特征子集以降低復雜度。本文還通過實例,闡明了共性特征提取方法、步驟與概念表達的方法,最后在給出的共性特征標示圖中驗證了該方法的效果和可行性。(4)顏色直方圖提取顏色特征顏色直方圖是在許多圖像檢索系統(tǒng)中被廣泛采用的顏色特征，具有特征提取和相似計算簡便的優(yōu)點。Swain 和Ballard 首先把直方圖應用于圖像檢索1,Smith J.R.提出了在HSV（Hue，Saturation，Value）空間的166 色量化方法2.但在對彩色圖像進行檢索時，顏色直方圖還存在兩個問題：（1）計算顏色直方圖需要將顏色空間劃分成若干個小的顏色區(qū)間，當區(qū)間的數(shù)目很大的時候不但會增加

20、計算負擔，也不利于在圖像庫中建立索引；（2）顏色直方圖丟棄了圖像的空間信息，對于兩幅顏色分布一樣，容并不一樣的圖像缺乏區(qū)分度。首先，根據(jù)HSV 顏色空間的色彩（Hue）、飽和度（Saturation）和亮度（Value）的值將圖像的顏色劃分為9 個區(qū)間，計算出每一種顏色區(qū)間所占的比重作為直方圖信息。然后計算出每一顏色區(qū)間所包含象素的空間分布信息。在圖像檢索時，綜合考慮顏色直方圖和空間信息的相似度。顏色直方圖特征，HSV 顏色空間有三個屬性：色彩（Hue），又稱為色調(diào)，是指色彩的相貌，與色彩明暗無關；飽和度（Saturation），是色彩純與不純的分別；亮度（Value），是指色彩的明暗程度，光

21、度的高低，越接近白色亮度越高，越接近灰色或黑色亮度越低。根據(jù)對HSV 空間的研究，人們能夠分辨出來的顏色大體上包括：紅、黃、綠、青、藍、紫，他們在色彩上的分布并非是均勻的，于是可以把色彩分為不等間隔的六份 Red，Yellow，Green，Cyan，Blue，Magenta，而亮度則可以劃分為不等間隔的黑灰白三份Black，Gray，White。這樣，HSV 顏色空間就劃分為9 個不同的空間區(qū)域。對于圖像中的一個象素，可以根據(jù)其Saturation 值的判定函數(shù)閾值將其劃分到相應顏色區(qū)間。3中采用了一個線性函數(shù) 1-0.8*Value 作為閾值，但是HSV 是一個不均勻的顏色空間，線性判別函數(shù)

22、不能準確地劃分彩色和非彩色的子空間?？梢圆捎靡环N更有效的非線性函數(shù)： 作為閾值。定義閾值：對于圖像中的任一象素X(H,S,V)，判斷其屬于哪一個顏色區(qū)間的公式由表1 所示：統(tǒng)計出每種顏色在圖像中占的比重就可以得到圖像的顏色直方圖信息。設圖像的顏色直方圖特征為His, 可由式（1）表示：顏色直方圖沒有表達出顏色空間分布的信息，也沒有反映出圖像中物體的形狀。因此，可以結合顏色區(qū)間的空間分布信息。顏色量化的結果也把圖像分為了9 個顏色區(qū)間，可以計算各個顏色區(qū)間所含象素的分布狀態(tài)作為其空間信息。設圖像顏色區(qū)間C 的空間信息為h(c)，它可以由公式（2）計算：顏色直方圖相似性計算：常用的基于顏色直方圖的

23、相似性度量方法包括直方圖相交、絕對值距離、歐式距離、加權距離等。本文采用直方圖相交作為圖像間的相似性度量標準。直方圖的相交是指兩個直方圖在每個顏色區(qū)間中共有的象素數(shù)量，度量方法如下：給定查詢圖像Q 和圖像庫中的圖像T，提取的圖像顏色直方圖特征分別為： q1 c , ,qNc c 和 t1 c , ,tNc c ，它們之間的相交距離由式（4）表示：空間信息的相似性計算：UFM4模糊匹配方法可以用來計算圖像的空間信息相似度，定義圖像顏色區(qū)間的平均距離為f d ，計算f d 的公式是（6）：其中，c N 為顏色區(qū)域數(shù)目。對于輸入圖像q P 和目標圖像t P ，它們的空間信息相似向量設為(Pq ,Pt

24、 ) M ,其計算公式：其中，i h 和j h 分別為輸入圖像和目標圖像的空間信息。設圖像顏色區(qū)間在圖像中的大小權重向量為w 。有多種選擇權重的方法，一種方法是每一個區(qū)間都賦予一樣的權重；另一種方法是以區(qū)間的面積在整個圖像中的比重作為其權重，即面積越大的區(qū)間越重要；還有一種方法是對位于圖像邊緣或者中心部分的區(qū)間賦予較大的權重。本文采用第二種方法。輸入圖像q P 和目標圖像t P 的空間信息相似度為k ，k 由公式（8）計算。綜合顏色直方圖和空間信息的圖像相似度：顏色直方圖雖然有計算簡單、對平移和旋轉不敏感的優(yōu)點，但無法捕捉顏色組成之間的空間關系，這無疑是不夠準確的。顏色直方圖量化本身就是一種粗

25、略的劃分方法，單獨采用基于顏色直方圖量化劃分的空間信息，檢索結果也不理想。本文采用綜合兩種特征的方法，對顏色直方圖和空間信息的相似度以不同的權重累加，檢索效果更加理想。設輸入圖像q P 和目標圖像t P 的相似度為 ， 由公式（9）計算：二、程序源代碼：2.1、邊界方向直方圖法提取形狀特征originailimage=imread('1.jpg'); image_I=rgb2gray(originailimage);% image_x,image_y=size(image_I);% if image_x<image_y% image_I=image_I' figu

26、re(1),subplot(121),title('灰度圖'), imshow(image_I); % 找到邊緣 I=image_I; image_bw=edge(I,'canny'); subplot(122),title('邊界檢測圖') imshow(image_bw); L,NUM=bwlabeln(image_bw,8); % 連通性 %將圖像邊緣量化為101份(1-100 101) 畫出邊緣直方圖，量化份數(shù)因情況不同而定 num_line=zeros(1,101); for i=1:NUM number=length(find(L=i

27、); for j=1:100 if number=j num_line(j)=num_line(j)+1; end end if number>100 num_line(101)=num_line(101)+1; end end EH_1=; for i=1:101 EH_1=EH_1 num_line(i); end x1=1:101; y1=EH_1; figure(2),title('plot曲線直方圖'), plot(x1,y1); %亦可將圖像量化為13份（1-10 11 12 13）畫出bar投影圖像 num=zeros(1,13); for i=1:NUM

28、number=length(find(L=i); for j=1:10 if number=j num(j)=num(j)+1; end end if number>10 & number<=20 num(11)=num(11)+1; end if number>20 & number<=100 num(12)=num(12)+1; end if number>100 num(13)=num(13)+1; end end EH=; for i=1:13 EH=EH num(i); end x=1:13; y=EH; figure(3),title(

29、'bar投影直方圖'), bar(x,y);2.2、基于PCA的圖像數(shù)據(jù)特征提取clear all% Define databasex=imread('3.jpg');figure(1)plot(x(:,1),x(:,2),'ko')xlabel('Valve Position (%)')ylabel('Fermenter Temperature (C)')% Calculate the covariance of dataset x (relationship between Valve and% Tempera

30、ture)covariance=cov(x);% Singular Value Decomposition x=T*E*P% E is the eignvalue of x; U is orthnormal; V=U'.%U E V=svd(covariance);%U=U' V=V'U E = eig(covariance);E = diag(E(2,2),E(1,1); % arrange the Eignvalue in a decending orderU = U(:,2),U(:,1); % arrange the Eignvector to make it

31、compatable with Eignvalue.% Mean center the raw dataValvebar=x(:,1)-mean(x(:,1);Tempbar=x(:,2)-mean(x(:,2);xbar=Valvebar,Tempbar;% Principal axis rotation of the covariance matrix.z=U'*xbar'z=z'%covz is the covariance matrix of principal components which is equal to E%matrix. The varianc

32、e of transformed variable (z1 & z2) will have variance%0.2283 and 0.0123 respectively.covz=U'*covariance*U;% The first column of U factor is - and + with nearly equal value which% means the first principal component is related to variability which both% measurements have difference. The seco

33、nd column is both positive which% means the 2nd PC is concentrate on the varaiability of the common between% the two.% scaling of PCS.for i=1:2Vs(:,i)=sqrt(E(i,i)*U(:,i);Ws(:,i)=U(:,i)/sqrt(E(i,i);end% after rescale: Vs'*Vs=E; V'*COV*Vs=E2; Ws'*Ws=inv(E); Ws'*COV*Ws=Iy=Ws'*xbar&#

34、39;%define T2T2=diag(y'*y);% find the control limit of T2p=2; %PC numbern=10; % sample numberT2Upper=p*(n-1)/(n-p)*finv(0.95,p,n-p);% draw the T2 control chartsfigure(2)plot(T2,'ko')hold onplot(T2Upper*ones(1,13),'k')xlabel('Sample Number')ylabel('T2 Score')hold o

35、ff% Draw the control ellipse% s12*s22/(s12*s22-s122)*(x1-x1mean)2/s12+(x2-x2mean)2/s22-2*s12(% x1-x1mean)(x2-x2mean)/s12/s22=T22(upper control limit)% S=covariance;% a=S(1)*S(4)/(S(1)*S(4)-S(2)2);% b=xbar(:,1).2/S(1)+xbar(:,2).2/S(4)-2*S(2)*xbar(:,1).*xbar(:,2)./S(1)/S(4);% error detection% two new

36、points is added to be detected which are 10.1,25.2;9.9;24;% plot two values with the original 10.Valve2=10.2,9.9;Temp2=25.1,25.5;figure(1)hold onplot(Valve2,Temp2,'k*')hold off% plot the control charts for valve and temperature separately.figure(3)hold onplot(1:10,Valve,'ko')plot(11:

37、12,Valve2,'k*')plot(0:13, mean(Valve)*ones(1,14),'k')plot(0:13, (mean(Valve)-std(Valve)*2)*ones(1,14),'k')plot(0:13, (mean(Valve)+std(Valve)*2)*ones(1,14),'k')xlabel('Sample Number')ylabel('Valve Position (%)')hold offfigure(4)hold onplot(1:10,Temp,

38、9;ko')plot(11:12,Temp2,'k*')plot(0:13, mean(Temp)*ones(1,14),'k')plot(0:13, (mean(Temp)+std(Temp)*2)*ones(1,14),'k')plot(0:13, (mean(Temp)-std(Temp)*2)*ones(1,14),'k')xlabel('Sample Number')ylabel('Temperature (C)')hold off% transform new observati

39、ons to PC axisy2=Ws'*(Valve2-mean(Valve);Temp2-mean(Temp);% plot the control charts for principal componentsfigure(5)hold onplot(1:10,y(1,:),'ko')plot(11:12,y2(1,:),'k*')plot(0:13, mean(y(1,:)*ones(1,14),'k')plot(0:13, -2*ones(1,14),'k')plot(0:13, 2*ones(1,14),

40、9;k')xlabel ('Sample Number')ylabel ('Principal Component 1')hold offfigure(6)hold onplot(1:10,y(2,:),'ko')plot(11:12,y2(2,:),'k*')plot(0:13,mean(y(2,:)*ones(1,14),'k')plot(0:13, -2*ones(1,14),'k')plot(0:13, 2*ones(1,14),'k')xlabel ('Sa

41、mple Number')ylabel ('Principal Component 2')hold off% plot the T2 chart for the two new pointsfigure(2)T22=diag(y2'*y2);hold onplot(11,12,T22,'k*')hold off% draw the control ellipsera = sqrt(T2Upper*E(1,1);rb = sqrt(T2Upper*E(2,2);V=sqrt(E(1,1)*U(:,1),sqrt(E(2,2)*-U(:,2);k1

42、= V(2,1)/V(1,1);k2 = V(2,2)/V(1,2);ang = atan(k1); %atan(k1);x0 = mean(x(:,1);y0 = mean(x(:,2);% i = 9.5:0.01:10.5;% s = 9.72:0.01:10.3;% j = k1*i + 42.6303;% l = k2*s + 19.3279;figure (7)h=ellipse(ra,rb,ang,x0,y0,'k',1000);hold onplot(x(:,1),x(:,2),'ko')plot(Valve2, Temp2,'k*

43、9;)% plot(i,j)% plot(s,l)xlabel ('Valve Position')ylabel ('Temperature')hold off2.3、基于Tamura紋理特征的紋理特征提取%計算Tamura紋理特征close allclear allclct0=cputime;I = imread('F-18.bmp');Nx,Ny = size(I);Ng=256;G=double(I);%計算粗糙度(coarseness)Sbest=zeros(Nx,Ny);E0h=zeros(Nx,Ny);E0v=zeros(Nx,Ny)

44、;E1h=zeros(Nx,Ny);E1v=zeros(Nx,Ny);E2h=zeros(Nx,Ny);E2v=zeros(Nx,Ny);E3h=zeros(Nx,Ny);E3v=zeros(Nx,Ny);E4h=zeros(Nx,Ny);E4v=zeros(Nx,Ny);E5h=zeros(Nx,Ny);E5v=zeros(Nx,Ny);flag=0;for i=1:Nx for j=2:Ny E0h(i,j)=G(i,j)-G(i,j-1); endendE0h=E0h/2;for i=1:Nx-1 for j=1:Ny E0v(i,j)=G(i,j)-G(i+1,j); endendE0

45、v=E0v/2;%圖片大小必須大于4*4才能計算E1h、E1vif (Nx<4|Ny<4) flag=1;endif(flag=0) for i=1:Nx-1 for j=3:Ny-1 E1h(i,j)=sum(sum(G(i:i+1,j:j+1)-sum(sum(G(i:i+1,j-2:j-1); end end for i=2:Nx-2 for j=2:Ny E1v(i,j)=sum(sum(G(i-1:i,j-1:j)-sum(sum(G(i+1:i+2,j-1:j); end end E1h=E1h/4; E1v=E1v/4;end%圖片大小必須大于8*8才能計算E2h、E

46、2vif (Nx<8|Ny<8) flag=1;endif(flag=0) for i=2:Nx-2 for j=5:Ny-3 E2h(i,j)=sum(sum(G(i-1:i+2,j:j+3)-sum(sum(G(i-1:i+2,j-4:j-1); end end for i=4:Nx-4 for j=3:Ny-1 E2v(i,j)=sum(sum(G(i-3:i,j-2:j+1)-sum(sum(G(i+1:i+4,j-2:j+1); end end E2h=E2h/16; E2v=E2v/16;end%圖片大小必須大于16*16才能計算E3h、E3vif (Nx<16|

47、Ny<16) flag=1endif(flag=0) for i=4:Nx-4 for j=9:Ny-7 E3h(i,j)=sum(sum(G(i-3:i+4,j:j+7)-sum(sum(G(i-3:i+4,j-8:j-1); end end for i=8:Nx-8 for j=5:Ny-3 E3v(i,j)=sum(sum(G(i-7:i,j-4:j+3)-sum(sum(G(i+1:i+8,j-4:j+3); end end E3h=E3h/64; E3v=E3v/64;end%圖片大小必須大于32*32才能計算E4h、E4v if (Nx<32|Ny<32) fla

48、g=1;endif(flag=0) for i=8:Nx-8 for j=17:Ny-15 E4h(i,j)=sum(sum(G(i-7:i+8,j:j+15)-sum(sum(G(i-7:i+8,j-16:j-1); end end for i=16:Nx-16 for j=9:Ny-7 E4v(i,j)=sum(sum(G(i-15:i,j-8:j+7)-sum(sum(G(i+1:i+16,j-8:j+7); end end E4h=E4h/256; E4v=E4v/256;end%圖片大小必須大于64*64才能計算E5h、E5vif (Nx<64|Ny<64) flag=1

49、;endif(flag=0) for i=16:Nx-16 for j=33:Ny-31 E5h(i,j)=sum(sum(G(i-15:i+16,j:j+31)-sum(sum(G(i-15:i+16,j-32:j-31); end end for i=32:Nx-32 for j=17:Ny-15 E5v(i,j)=sum(sum(G(i-31:i,j-16:j+15)-sum(sum(G(i+1:i+32,j-16:j+15); end end E5h=E5h/1024; E5v=E5v/1024;endfor i=1:Nx for j=1:Ny maxv,index=max(E0h(i

50、,j),E0v(i,j),E1h(i,j),E1v(i,j),E2h(i,j),E2v(i,j),E3h(i,j),E3v(i,j),E4h(i,j),E4v(i,j),E5h(i,j),E5v(i,j); k=floor(index+1)/2); Sbest(i,j)=2.k; endendFcoarseness=sum(sum(Sbest)/(Nx*Ny);%計算對比度counts,graylevels=imhist(I);PI=counts/(Nx*Ny);averagevalue=sum(graylevels.*PI);u4=sum(graylevels-repmat(averagev

51、alue,256,1).4.*PI);standarddeviation=sum(graylevels-repmat(averagevalue,256,1).2.*PI);alpha4=u4/standarddeviation2;Fcontrast=sqrt(standarddeviation)/alpha4.(1/4);%計算方向度PrewittH=-1 0 1;-1 0 1;-1 0 1;PrewittV=1 1 1;0 0 0;-1 -1 -1;%計算橫向梯度deltaH=zeros(Nx,Ny);for i=2:Nx-1 for j=2:Ny-1 deltaH(i,j)=sum(sum(G(i-1:i+1,j-1:j+1).*PrewittH); endendfor j=2:Ny-1 deltaH(1,j)=G(1,j+1)-G(1,j); deltaH(Nx,j)=G(Nx,j+1)-G(Nx,j); endfor i=1:Nx deltaH(i,1)