PeakFinder

Функции
void	PeakFinder::findPeaks (std::vector< float > x0, std::vector< int > &peakInds, bool includeEndpoints=true, float extrema=1)

Подробное описание

Функции

findPeaks()

void PeakFinder::findPeaks	(	std::vector< float >	x0,
		std::vector< int > &	peakInds,
		bool	includeEndpoints = true,
		float	extrema = 1
	)

Аргументы

x0	входной массив
extrema	1 - если нужен максимум, -1 - если нужен минимум
includeEndpoints	true - если границы интервала входного массива x0 должны быть включены в поиск пика, false - в противном случае

Возвращает

peakInds: индексы пиков, найденные в x0

См. определение в файле peak_finder.cpp строка 66

{
    int minIdx = distance(x0.begin(), min_element(x0.begin(), x0.end()));
    int maxIdx = distance(x0.begin(), max_element(x0.begin(), x0.end()));
 
    float sel = (x0[maxIdx]-x0[minIdx])/4.0;
    int len0 = x0.size();
 
    scalarProduct(extrema, x0, x0);
 
    std::vector<float> dx;
    diff(x0, dx);
    replace(dx.begin(), dx.end(), 0.0f, -PeakFinder::EPS);
    std::vector<float> dx0(dx.begin(), dx.end()-1);
    std::vector<float> dx0_1(dx.begin()+1, dx.end());
    std::vector<float> dx0_2;
 
    vectorElementsProduct(dx0, dx0_1, dx0_2);
 
    std::vector<int> ind;
    findIndicesLessThan(dx0_2, 0, ind); // Find where the derivative changes sign   
    std::vector<float> x;
    float leftMin;
    int minMagIdx;
    float minMag;
    
    if(includeEndpoints)
    {
        //x = [x0(1);x0(ind);x0(end)];  
        selectElementsFromIndices(x0, ind, x);      
        x.insert(x.begin(), x0[0]);
        x.insert(x.end(), x0[x0.size()-1]);
        //ind = [1;ind;len0];
        ind.insert(ind.begin(), 1);
        ind.insert(ind.end(), len0);
        minMagIdx = distance(x.begin(), std::min_element(x.begin(), x.end()));
        minMag = x[minMagIdx];      
        //std::cout<<"Hola"<<std::endl;
        leftMin = minMag;
    }
    else
    {
        selectElementsFromIndices(x0, ind, x);
        if(x.size()>2)
        {
            minMagIdx = distance(x.begin(), std::min_element(x.begin(), x.end()));      
            minMag = x[minMagIdx];              
            leftMin = x[0]<x0[0]?x[0]:x0[0];
        }
    }
 
    int len = x.size();
 
    if(len>2)
    {
        float tempMag = minMag;
        bool foundPeak = false;
        int ii;
 
        if(includeEndpoints)
        {
            // Deal with first point a little differently since tacked it on
            // Calculate the sign of the derivative since we tacked the first
            //  point on it does not neccessarily alternate like the rest.
            std::vector<float> xSub0(x.begin(), x.begin()+3);//tener cuidado subvector
            std::vector<float> xDiff;//tener cuidado subvector
            diff(xSub0, xDiff);
 
            std::vector<int> signDx;
            signVector(xDiff, signDx);
 
            if (signDx[0] <= 0) // The first point is larger or equal to the second
            {
                if (signDx[0] == signDx[1]) // Want alternating signs
                {
                    x.erase(x.begin()+1);
                    ind.erase(ind.begin()+1);
                    len = len-1;
                }
            }
            else // First point is smaller than the second
            {
                if (signDx[0] == signDx[1]) // Want alternating signs
                {
                    x.erase(x.begin());
                    ind.erase(ind.begin());
                    len = len-1;
                }
            }
        }
 
        //Skip the first point if it is smaller so we always start on maxima
        if ( x[0] >= x[1] )
            ii = 0;
        else
            ii = 1;
 
        //Preallocate max number of maxima
        float maxPeaks = ceil((float)len/2.0);
        std::vector<int> peakLoc(maxPeaks,0);
        std::vector<float> peakMag(maxPeaks,0.0);
        int cInd = 1;
        int tempLoc = 0;        
    
        while(ii < len)
        {
            ii = ii+1;//This is a peak
            //Reset peak finding if we had a peak and the next peak is bigger
            //than the last or the left min was small enough to reset.
            if(foundPeak)
            {
                tempMag = minMag;
                foundPeak = false;
            }
        
            //Found new peak that was lager than temp mag and selectivity larger
            //than the minimum to its left.
        
            if( x[ii-1] > tempMag && x[ii-1] > leftMin + sel )
            {
                tempLoc = ii-1;
                tempMag = x[ii-1];
            }
 
            //Make sure we don't iterate past the length of our vector
            if(ii == len)
                break; //We assign the last point differently out of the loop
 
            ii = ii+1; // Move onto the valley
            
            //Come down at least sel from peak
            if(!foundPeak && tempMag > sel + x[ii-1])
            {               
                foundPeak = true; //We have found a peak
                leftMin = x[ii-1];
                peakLoc[cInd-1] = tempLoc; // Add peak to index
                peakMag[cInd-1] = tempMag;
                cInd = cInd+1;
            }
            else if(x[ii-1] < leftMin) // New left minima
                leftMin = x[ii-1];
            
        }
 
        // Check end point
        if(includeEndpoints)
        {        
            if ( x[x.size()-1] > tempMag && x[x.size()-1] > leftMin + sel )
            {
                peakLoc[cInd-1] = len-1;
                peakMag[cInd-1] = x[x.size()-1];
                cInd = cInd + 1;
            }
            else if( !foundPeak && tempMag > minMag )// Check if we still need to add the last point
            {
                peakLoc[cInd-1] = tempLoc;
                peakMag[cInd-1] = tempMag;
                cInd = cInd + 1;
            }
        }
        else if(!foundPeak)
        {
            float minAux = x0[x0.size()-1]<x[x.size()-1]?x0[x0.size()-1]:x[x.size()-1];
            if ( x[x.size()-1] > tempMag && x[x.size()-1] > leftMin + sel )
            {
                peakLoc[cInd-1] = len-1;
                peakMag[cInd-1] = x[x.size()-1];
                cInd = cInd + 1;
            }
            else if( !tempMag >  minAux + sel)// Check if we still need to add the last point
            {
                peakLoc[cInd-1] = tempLoc;
                peakMag[cInd-1] = tempMag;
                cInd = cInd + 1;
            }
        }
 
        //Create output
        if( cInd > 0 )
        {           
            std::vector<int> peakLocTmp(peakLoc.begin(), peakLoc.begin()+cInd-1);
            selectElementsFromIndices(ind, peakLocTmp, peakInds);           
        }       
 
    }
    //else
    //{
        //input signal length <= 2
    //}
}