1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
|
// ************************************************************************************************
//
// BornAgain: simulate and fit reflection and scattering
//
//! @file Device/Analyze/Peaks.cpp
//! @brief Implements function FindPeaks namespace Analyze.
//!
//! @homepage http://www.bornagainproject.org
//! @license GNU General Public License v3 or higher (see COPYING)
//! @copyright Forschungszentrum Jülich GmbH 2018
//! @authors Scientific Computing Group at MLZ (see CITATION, AUTHORS)
//
// ************************************************************************************************
#include "Device/Analyze/Peaks.h"
#include "Base/Axis/Scale.h"
#include "Device/Data/DataUtil.h"
#include "Device/Data/Datafield.h"
#include <tspectrum.h> // third-party code, extracted from CERN ROOT (class TSpectrum2)
std::vector<std::pair<double, double>>
Analyze::FindPeaks(const Datafield& data, double sigma, const std::string& option, double threshold)
{
const double2d_t arr = data.values2D();
tspectrum::Spectrum2D spec;
const auto peaks = spec.find_peaks(arr, sigma, option, threshold);
// coordinates of peaks in histogram axes units
std::vector<std::pair<double, double>> result;
for (const auto& p : peaks) {
const double row_value = p.first;
const double col_value = p.second;
const auto xaxis_index = static_cast<size_t>(col_value);
const size_t yaxis_index = data.yAxis().size() - 1 - static_cast<size_t>(row_value);
const Bin1D xbin = data.xAxis().bin(xaxis_index);
const Bin1D ybin = data.yAxis().bin(yaxis_index);
const double dx = col_value - static_cast<size_t>(col_value);
const double dy = -1.0 * (row_value - static_cast<size_t>(row_value));
const double x = xbin.center() + xbin.binSize() * dx;
const double y = ybin.center() + ybin.binSize() * dy;
result.emplace_back(x, y);
}
return result;
}
|