Cepstral Pitch Tracker

//    std::cerr << "findPeakLocation: size " << size << ", peakIndex " << peakIndex << std::endl;

//        std::cerr << "returning " << peakIndex << ", data too short" << std::endl;

    double maxval = 0.0;

    double location = peakIndex;

    const int divisions = 10;

    double y[4];

    y[0] = data[peakIndex-1];

    y[1] = data[peakIndex];

    y[2] = data[peakIndex+1];

    if (peakIndex < size - 2) {

        y[3] = data[peakIndex+2];

    } else {

        y[3] = y[2];

    }

//    std::cerr << "a y: " << y[0] << " " << y[1] << " " << y[2] << " " << y[3] << std::endl;

    for (int i = 0; i < divisions; ++i) {

        double probe = double(i) / double(divisions);

        double value = cubicInterpolate(y, probe);

//        std::cerr << "probe = " << probe << ", value = " << value << " for location " << peakIndex + probe << std::endl;

        if (value > maxval) {

            maxval = value; 

            location = peakIndex + probe;

        }

    y[3] = y[2];

    y[2] = y[1];

    y[1] = y[0];

    if (peakIndex > 1) {

        y[0] = data[peakIndex-2];

    } else {

        y[0] = y[1];

    }

//    std::cerr << "b y: " << y[0] << " " << y[1] << " " << y[2] << " " << y[3] << std::endl;

    for (int i = 0; i < divisions; ++i) {

        double probe = double(i) / double(divisions);

        double value = cubicInterpolate(y, probe);

//        std::cerr << "probe = " << probe << ", value = " << value << " for location " << peakIndex - 1 + probe << std::endl;

        if (value > maxval) {

            maxval = value; 

            location = peakIndex - 1 + probe;

        }

    }

//    std::cerr << "returning " << location << std::endl;

    return location;