投影曲線的波峰、波谷查找

一、前言

在圖像分析里，投影曲線是我們經常要用到的一個圖像特征，通過投影曲線我們可以看到在某一個方向上，圖像灰度變化的規律，這在圖像分割，文字提取方面應用比較廣。一個投影曲線，它的關鍵信息就在於波峰與波谷，所以我們面臨的第一個問題就是找到波峰與波谷。

第一次涉及到求波峰與波谷時，很多人都不以為意，覺得波谷波峰還不容易，無非是一些曲線變化為零的點，從離散的角度來說，也就是:

波峰：F(x)>F(x−1)且F(x)>F(x+1)

波谷：F(x)<F(x−1)且F(x)<F(x+1)

這么簡單嗎？顯示不是，你首先就會遇到這樣的曲線圖，然后圖上的波峰點並不滿足上面的條件。

看到這種情況，你也許會考慮在上面的等式中把>和<改為≥和≤。

波峰：F(x)≥F(x−1)&&F(x)>F(x+1)  或者 F(x)>F(x−1)&&F(x)≥F(x+1)

波谷：F(x)≤F(x−1)&&F(x)<F(x+1)  或者 F(x)<F(x−1)&&F(x)≤F(x+1)

這次是否就這樣簡單，答案顯示不是，下面的這個圖就會讓你對一些非峰值點作出錯誤的判斷。

上面這幅圖真正的峰值只有一個，其他平台上的點，你如果按上面修改的公式，就會被錯誤的當成波峰點。

下面讓我們看一下，到底如何能求得准確的曲線波峰與波谷。

二、波峰波谷算法

投影曲線實際上是一個一維的向量：

V=[v1,v2,…,vn]

 

其中vi,i∈[1,2,…,N]，代表圖像在第i行或列上的灰度累積。當然不僅僅是投影曲線，T也可以是某一事件中變量的觀測值，我們需要研究這個變量的變化規律。

下面給出波峰與波谷的算法：

1，假投影曲線可以表示為V=[v1,v2,…,vn]。

2，計算V的一階差分向量DiffV:

Diffv(i)=V(i+1)−V(i),其中i∈1,2,…,N−1

 

3,對差分向量進行取符號函數運算，Trend=sign(Diffv),即遍歷Diffv，若Diffv(i)大於0，則取1；如果小於0，則取-1，否則則值為0。

sign(x)=⎧⎩⎨⎪⎪10−1if  x>0if  x=0if  x<0

 

4，從尾部遍歷Trend向量,進行如下操作：

if Trend(i)=0且Trend(i+1)≥0，則Trend(i)=1if Trend(i)=0且Trend(i+1)<0，則Trend(i)=−1

 

5，對Trend向量進行一階差分運算，如同步驟2，得到R=diff(Trend)。

6，遍歷得到的差分向量R，如果R(i)=−2，則i+1為投影向量V的一個峰值位，對應的峰值為V(i+1)；如果R(i)=2，則i+1為投影向量V的一個波谷位，對應的波谷為V(i+1)。

下面我們來結合一個實際的向量值，給中中間結合的計算。

1，V=[−5,10,10,14,14,8,8,6,6,−3,2,2,2,2,−3]。

它的曲線圖像如下把示，圖中紅色圈標出了曲線的峰值，而綠字圈標出了圖像的波谷位置。

2，計算V的一階差分，我們得到Diff(V)=[15,0,4,0,,−6,0,−2,0,−9,5,0,0,0,−5]。

3，對Diffv進行取符號運算，得到向量Trend=[1,0,1,0,−1,0,−1,0,−1,1,0,0,0,1]。

4，對Trend作一次遍歷，如步驟4。Trend=[1,1,1,−1,−1,−1,−1,−1,−1,1,−1,−1,−1]。

5，對Trend做一階差分，得到向量R=Diff(Trend)=[0,0,−2,,0,0,0,0,0,2,−2,0,0,0]。

6，遍歷向量R,我們就得到了兩個峰值點和一個波谷點。

三、算法原理

其實上述算法的核心思路非常簡單，曲線的峰值點，滿足一階導數為0，並且滿足二階導數為負；而波谷點，則滿足一階導數為0，二階導數為正。

在上面的算法里面，我們首先計算了一階的導數Diffv，然后我們將其符號化，是因為我們並不關心一階導數的大小。

然后我們去看那些一階層數為0的地方，我們發現，那些平台上的點，有些並不是波峰與波谷，然后很多處在上坡與下坡的路上，所以我們將它們的一階導數設為與它們所在的坡面梯度方向相同。

最后我們再來計算二階導數時，就會發現只要為2或者-2,所以曲線斜在這個點發生了變化，由正變負或由負變正。找到這些點，也就找到了原曲線中的波峰或波谷點。

四、實現

下面給出這個算法的C++實現，findPeaks是查找波峰函數，而查找波谷函數則類似，這里沒有寫在一個函數內。函數接受一個Vecotr<int>的向量，輸出為一個vector<int>的位置向量。

 1 void findPeak(const vector<int>& v, vector<int>& peakPositions)  2 {  3 vector<int> diff_v (v.size() - 1, 0);  4 // 計算V的一階差分和符號函數trend  5 for (vector<int>::size_type i = 0; i !=diff_v.size(); i++)  6  {  7 if (v[i + 1] - v[i]>0)  8 diff_v[i] = 1;  9 else if (v[i + 1] - v[i] < 0) 10 diff_v[i] = -1; 11 else 12 diff_v[i] = 0; 13  } 14 // 對Trend作了一個遍歷 15 for (int i = diff_v.size() - 1; i >= 0; i--) 16  { 17 if (diff_v[i] == 0 && i ==diff_v.size()-1) 18  { 19 diff_v[i] = 1; 20  } 21 else if (diff_v[i] == 0) 22  { 23 if (diff_v[i + 1] >= 0) 24 diff_v[i] = 1; 25 else 26 diff_v[i] = -1; 27  } 28  } 29 30 for (vector<int>::size_type i = 0; i != diff_v.size() - 1; i++) 31  { 32 if (diff_v[i + 1] - diff_v[i] == -2) 33 peakPositions.push_back(i + 1); 34  } 35 }