Skip to content
Snippets Groups Projects
Commit 405d13cd authored by Vickie Lynch's avatar Vickie Lynch
Browse files

Refs #7449 struct for sorting peaks

parent 808795a5
No related branches found
No related tags found
No related merge requests found
...@@ -54,6 +54,19 @@ namespace Algorithms ...@@ -54,6 +54,19 @@ namespace Algorithms
/// Implement abstract Algorithm methods /// Implement abstract Algorithm methods
void exec(); void exec();
double moment(MantidVec& X, size_t n, double mean, int k); double moment(MantidVec& X, size_t n, double mean, int k);
struct cont_peak
{
size_t start;
size_t stop;
double maxY;
};
struct by_len
{
bool operator()(cont_peak const &a, cont_peak const &b)
{
return a.maxY > b.maxY;
}
};
}; };
......
...@@ -147,54 +147,58 @@ namespace Algorithms ...@@ -147,54 +147,58 @@ namespace Algorithms
} }
while (std::abs(Ymean-Yvariance) > k * Ysigma); while (std::abs(Ymean-Yvariance) > k * Ysigma);
// remove single outliers if(n > 5)
{
if (mask[1] == mask[2] && mask[2] == mask[3]) // remove single outliers
mask[0] = mask[1]; if (mask[1] == mask[2] && mask[2] == mask[3])
if (mask[0] == mask[2] && mask[2] == mask[3]) mask[0] = mask[1];
mask[1] = mask[2]; if (mask[0] == mask[2] && mask[2] == mask[3])
for (size_t l = 2; l < n-2; ++l) mask[1] = mask[2];
{ for (size_t l = 2; l < n-3; ++l)
if (mask[l-1] == mask[l+1] && (mask[l-1] == mask[l-2] || mask[l+1] == mask[l+2]))
{
mask[l] = mask[l+1];
}
}
if (mask[n-2] == mask[n-3] && mask[n-3] == mask[n-4])
mask[n-1] = mask[n-2];
if (mask[n-1] == mask[n-3] && mask[n-3] == mask[n-4])
mask[n-2] = mask[n-1];
// mask regions not connected to largest region
vector<size_t> cont_start;
vector<size_t> cont_stop;
for (size_t l = 1; l < n-1; ++l)
{
if (mask[l] != mask[l-1] && mask[l] == 1)
{ {
cont_start.push_back(l); if (mask[l-1] == mask[l+1] && (mask[l-1] == mask[l-2] || mask[l+1] == mask[l+2]))
{
mask[l] = mask[l+1];
}
} }
if (mask[l] != mask[l-1] && mask[l] == 0) if (mask[n-2] == mask[n-3] && mask[n-3] == mask[n-4])
mask[n-1] = mask[n-2];
if (mask[n-1] == mask[n-3] && mask[n-3] == mask[n-4])
mask[n-2] = mask[n-1];
// mask regions not connected to largest region
// for loop can start > 1 for multiple peaks
vector<cont_peak> peaks;
for (size_t l = 1; l < n; ++l)
{ {
cont_stop.push_back(l-1); if (mask[l] != mask[l-1] && mask[l] == 1)
{
peaks.push_back(cont_peak());
peaks[peaks.size()-1].start = l;
}
if (peaks.size() > 0)
{
size_t ipeak = peaks.size()-1;
if (mask[l] != mask[l-1] && mask[l] == 0)
{
peaks[ipeak].stop = l-1;
}
if (inpY[l] > peaks[ipeak].maxY) peaks[ipeak].maxY = inpY[l];
}
} }
} if(peaks.size()> 0)
if(cont_start.size() > cont_stop.size()) cont_stop.push_back(n-1);
vector<size_t> cont_len;
for (size_t l = 0; l < cont_start.size(); ++l)
{
cont_len.push_back(cont_stop[l] - cont_start[l] + 1);
}
std::vector<size_t>::iterator ic = std::max_element(cont_len.begin(), cont_len.end());
const size_t c = ic - cont_len.begin();
for (size_t l = 0; l < cont_len.size(); ++l)
{
if (l != c) for (size_t j = cont_start[l]; j <= cont_stop[l]; ++j)
{ {
mask[j] = 0; if(peaks[peaks.size()-1].stop == 0) peaks[peaks.size()-1].stop = n-1;
std::sort(peaks.begin(), peaks.end(), by_len());
for (size_t l = 1; l < peaks.size(); ++l)
{
for (size_t j = peaks[l].start; j <= peaks[l].stop; ++j)
{
mask[j] = 0;
}
}
} }
} }
// save output of mask * Y // save output of mask * Y
vecX = inpX; vecX = inpX;
std::transform(mask.begin(), mask.end(), inpY.begin(), vecY.begin(), std::multiplies<double>()); std::transform(mask.begin(), mask.end(), inpY.begin(), vecY.begin(), std::multiplies<double>());
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment