Revision 40:8f56ef28b0b1 PeakInterpolator.cpp
| PeakInterpolator.cpp | ||
|---|---|---|
| 24 | 24 |
|
| 25 | 25 |
#include "PeakInterpolator.h" |
| 26 | 26 |
|
| 27 |
#include <iostream> |
|
| 28 |
|
|
| 27 | 29 |
static double cubicInterpolate(const double y[4], double x) |
| 28 | 30 |
{
|
| 29 | 31 |
double a0 = y[3] - y[2] - y[0] + y[1]; |
| ... | ... | |
| 40 | 42 |
double |
| 41 | 43 |
PeakInterpolator::findPeakLocation(const double *data, int size, int peakIndex) |
| 42 | 44 |
{
|
| 43 |
if (peakIndex < 2 || peakIndex > size - 3) {
|
|
| 45 |
std::cerr << "findPeakLocation: size " << size << ", peakIndex " << peakIndex << std::endl; |
|
| 46 |
|
|
| 47 |
if (peakIndex < 1 || peakIndex > size - 2) {
|
|
| 48 |
std::cerr << "returning " << peakIndex << ", data too short" << std::endl; |
|
| 44 | 49 |
return peakIndex; |
| 45 | 50 |
} |
| 46 | 51 |
|
| ... | ... | |
| 53 | 58 |
y[0] = data[peakIndex-1]; |
| 54 | 59 |
y[1] = data[peakIndex]; |
| 55 | 60 |
y[2] = data[peakIndex+1]; |
| 56 |
y[3] = data[peakIndex+2]; |
|
| 61 |
if (peakIndex < size - 2) {
|
|
| 62 |
y[3] = data[peakIndex+2]; |
|
| 63 |
} else {
|
|
| 64 |
y[3] = y[2]; |
|
| 65 |
} |
|
| 66 |
std::cerr << "a y: " << y[0] << " " << y[1] << " " << y[2] << " " << y[3] << std::endl; |
|
| 57 | 67 |
for (int i = 0; i < divisions; ++i) {
|
| 58 | 68 |
double probe = double(i) / double(divisions); |
| 59 | 69 |
double value = cubicInterpolate(y, probe); |
| 70 |
std::cerr << "probe = " << probe << ", value = " << value << " for location " << peakIndex + probe << std::endl; |
|
| 60 | 71 |
if (value > maxval) {
|
| 61 | 72 |
maxval = value; |
| 62 | 73 |
location = peakIndex + probe; |
| ... | ... | |
| 66 | 77 |
y[3] = y[2]; |
| 67 | 78 |
y[2] = y[1]; |
| 68 | 79 |
y[1] = y[0]; |
| 69 |
y[0] = data[peakIndex-2]; |
|
| 80 |
if (peakIndex > 1) {
|
|
| 81 |
y[0] = data[peakIndex-2]; |
|
| 82 |
} else {
|
|
| 83 |
y[0] = y[1]; |
|
| 84 |
} |
|
| 85 |
std::cerr << "b y: " << y[0] << " " << y[1] << " " << y[2] << " " << y[3] << std::endl; |
|
| 70 | 86 |
for (int i = 0; i < divisions; ++i) {
|
| 71 | 87 |
double probe = double(i) / double(divisions); |
| 72 | 88 |
double value = cubicInterpolate(y, probe); |
| 89 |
std::cerr << "probe = " << probe << ", value = " << value << " for location " << peakIndex - 1 + probe << std::endl; |
|
| 73 | 90 |
if (value > maxval) {
|
| 74 | 91 |
maxval = value; |
| 75 | 92 |
location = peakIndex - 1 + probe; |
| 76 | 93 |
} |
| 77 | 94 |
} |
| 78 | 95 |
|
| 96 |
std::cerr << "returning " << location << std::endl; |
|
| 97 |
|
|
| 79 | 98 |
return location; |
| 80 | 99 |
} |
| 81 | 100 |
|
Also available in: Unified diff