Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
  • Loading branch information
FloCiaglia committed Nov 7, 2019
2 parents 0e569de + 0f15810 commit 8a6458d
Show file tree
Hide file tree
Showing 2 changed files with 94 additions and 21 deletions.
80 changes: 71 additions & 9 deletions src/GaussianFitter.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -614,19 +614,81 @@ int GaussianFitter::guess_peaks(std::vector<Peak*>* results,
//are pointing to space used in LidarVolume
results->clear();

//Sign of gradient:
// = 1 for increasing
// = 0 for level AND PREVIOUSLY INCREASING (so potential wide peak)
// = -1 for decreasing OR level, but previously decreasing
//A sharp peak is identified by grad=1 -> grad=-1
//A wide peak is identified by grad=0 -> grad=-1
std::vector<float> peak_guesses_loc; //Store peaks x-values here
std::vector<float> peak_guesses_amp;
int wideStart = -1; //The start of any current wide peak
int prev_grad = -1;
int grad = -1;
for(int i = 0; i<(int)ampData.size()-1; i++){
int firstDiffs[(int)ampData.size()];
int secondDiffs[(int)ampData.size()];
firstDiffs[1] = ampData[1] - ampData[0];
for(int i = 2; i<(int)ampData.size(); i++){
firstDiffs[i] = ampData[i] - ampData[i-1];
secondDiffs[i] = firstDiffs[i] - firstDiffs[i-1];
}
for(int i = 2; i<(int)ampData.size(); i++){
if(grad == 0){
if(firstDiffs[i] > 0){
// flat to positive (record peak if not a trough)
if(ampData[i]>noise_level && prev_grad == 1){
//handle flat peaks
int j = i-2;
for( ; ampData[j] == ampData[i-1] ; j--);
// there is no body to the above for loop -- just looking
// for the stopping point
float lenOfFlat = i-1 -j -1.;
//this is truncating, we need
//to fix this
float loc = i-1 - lenOfFlat/2.;
//record the peak
peak_guesses_loc.push_back(loc);
peak_guesses_amp.push_back(ampData[i-1]);
}
grad = 1;
}else if(firstDiffs[i] < 0 ){
// flat to negative (record peak)
if(ampData[i]>noise_level){
//handle flat peaks
int j = i-2;
for( ; ampData[j] == ampData[i-1] ; j--);
// there is no body to the above for loop -- just looking
// for the stopping point
float lenOfFlat = i-1 -j -1.;
//this is truncating, we need
//to fix this
float loc = i-1 - lenOfFlat/2.;
//record the peak
peak_guesses_loc.push_back(loc);
peak_guesses_amp.push_back(ampData[i-1]);
}
grad = -1;
}

}else if(grad == 1){
prev_grad = 1;
if(firstDiffs[i] == 0 ){
grad = 0;
}else if(firstDiffs[i] < 0){
// positive to negative (record peak)
if(ampData[i]>noise_level){
//record the peak
peak_guesses_loc.push_back(float(i-1));
peak_guesses_amp.push_back(ampData[i-1]);
}
grad = -1;
}

}else{ // gradient is -1
prev_grad = -1;
if(firstDiffs[i] == 0 ){
// down to flat
grad = 0;
}else if(firstDiffs[i] > 0){
// down to up
grad = 1;
}
}
}
/*for(int i = 0; i<(int)ampData.size()-1; i++){
// were we sloping up before?
if (grad == 1){
// sloping down
Expand Down Expand Up @@ -659,7 +721,7 @@ int GaussianFitter::guess_peaks(std::vector<Peak*>* results,
grad = 1;
}
}
}
}*/

//Figure out the size of ampData
size_t n = ampData.size();
Expand Down
35 changes: 23 additions & 12 deletions src/GaussianFitter_unittests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -72,15 +72,20 @@ TEST_F(GaussianFitterTest, NayaniClipped1_guess){
std::vector<Peak*> peaks;
int count = fitter.guess_peaks(&peaks, ampData, idxData);

EXPECT_EQ(2,peaks.size());
ASSERT_EQ(3,peaks.size());
EXPECT_EQ(200,peaks.at(0)->amp);
EXPECT_EQ(34 ,peaks.at(1)->amp);
EXPECT_NEAR(18.5, peaks.at(0)->location, .25);
EXPECT_EQ(38, peaks.at(1)->location);
EXPECT_NEAR(10.5, peaks.at(0)->fwhm, 1);
EXPECT_NEAR(6, peaks.at(1)->fwhm, 1);

EXPECT_EQ(2, count);
EXPECT_EQ(43 ,peaks.at(1)->amp);
EXPECT_NEAR(30.5, peaks.at(1)->location,.25);
EXPECT_NEAR(4, peaks.at(1)->fwhm, 1);

EXPECT_EQ(34 ,peaks.at(2)->amp);
EXPECT_EQ(38, peaks.at(2)->location);
EXPECT_NEAR(6, peaks.at(2)->fwhm, 1);

EXPECT_EQ(3, count);

}

Expand All @@ -102,7 +107,7 @@ TEST_F(GaussianFitterTest, NayaniClipped2_guess){

ASSERT_EQ(2,peaks.size());
EXPECT_EQ(235,peaks.at(0)->amp);
EXPECT_EQ(18, peaks.at(0)->location);
EXPECT_NEAR(17.5, peaks.at(0)->location,1);
EXPECT_NEAR(7.8, peaks.at(0)->fwhm, 1);

EXPECT_EQ(29, peaks.at(1)->location);
Expand Down Expand Up @@ -186,12 +191,15 @@ TEST_F(GaussianFitterTest, NayaniClipped4_guess){
std::vector<Peak*> peaks;
int count = fitter.guess_peaks(&peaks, ampData, idxData);

EXPECT_EQ(1,peaks.size());
ASSERT_EQ(2,peaks.size());
EXPECT_EQ(240,peaks.at(0)->amp);
EXPECT_EQ(20, peaks.at(0)->location);
EXPECT_NEAR(7.4, peaks.at(0)->fwhm, 1);

EXPECT_EQ(1, count);
EXPECT_EQ(16,peaks.at(1)->amp);
EXPECT_NEAR(29.5, peaks.at(1)->location,.25);
EXPECT_NEAR(8, peaks.at(1)->fwhm, 1);
EXPECT_EQ(2, count);

}

Expand Down Expand Up @@ -276,15 +284,18 @@ TEST_F(GaussianFitterTest, NayaniClipped7_guess){
fitter.noise_level = 9;
int count = fitter.guess_peaks(&peaks, ampData, idxData);

ASSERT_EQ(2,peaks.size());
ASSERT_EQ(3,peaks.size());
EXPECT_EQ(98,peaks.at(0)->amp);
EXPECT_EQ(168,peaks.at(1)->amp);
EXPECT_EQ(18, peaks.at(0)->location);
EXPECT_EQ(31, peaks.at(1)->location);
EXPECT_NEAR(5.2, peaks.at(0)->fwhm, 1);
EXPECT_NEAR(10, peaks.at(1)->fwhm, 1);
EXPECT_EQ(13,peaks.at(2)->amp);
EXPECT_NEAR(40.5, peaks.at(2)->location,.25);
EXPECT_NEAR(11, peaks.at(2)->fwhm, 1);

EXPECT_EQ(2, count);
EXPECT_EQ(3, count);

}

Expand Down Expand Up @@ -1681,10 +1692,10 @@ TEST_F(GaussianFitterTest, problem_waveform_5_find){
fitter.smoothing_expt(&ampData);
int count = fitter.find_peaks(&peaks,ampData,idxData, 200);

EXPECT_EQ(2,peaks.size());
ASSERT_EQ(2,peaks.size());
EXPECT_NEAR(19.5, peaks.at(0)->location,.25);
EXPECT_EQ(186, peaks.at(0)->amp);
EXPECT_EQ(14, peaks.at(1)->amp);
EXPECT_NEAR(19.5, peaks.at(0)->location,.25);
EXPECT_EQ(30, peaks.at(1)->location);
EXPECT_NEAR(5.7, peaks.at(0)->fwhm, 1);
EXPECT_NEAR(8, peaks.at(1)->fwhm, 1);
Expand Down

0 comments on commit 8a6458d

Please sign in to comment.