22namespace CurveFitting {
50 const std::string
name()
const override {
return "FitPowderDiffPeaks"; }
52 const std::string
summary()
const override {
return "Fit peaks in powder diffraction pattern. "; }
55 int version()
const override {
return 1; }
56 const std::vector<std::string>
seeAlso()
const override {
return {
"LeBailFit"}; }
59 const std::string
category()
const override {
return "Diffraction\\Fitting"; }
68 void processInputProperties();
75 std::map<std::string, double> parammap,
76 std::map<std::string, std::string> bk2bk2braggmap,
bool &good,
77 std::vector<int> &hkl,
double &d_h);
80 bool getHKLFromMap(std::map<std::string, int> intmap, std::vector<int> &hkl);
87 std::vector<std::map<std::string, double>> ¶mmaps,
88 std::vector<std::map<std::string, int>> &hklmaps);
91 void fitPeaksWithGoodStartingValues();
94 void fitPeaksRobust();
98 double leftdev,
double rightdev,
size_t m_wsIndex,
double &chi2);
105 double peakrightbound,
const std::map<std::string, double> &rightpeakparammap,
110 const std::vector<std::string> ¶mtodomc);
118 double rightbound,
double &chi2,
bool &annhilatedpeak);
122 double dampingfactor);
125 bool fitOverlappedPeaks(std::vector<Functions::BackToBackExponential_sptr> peaks,
131 std::vector<Functions::BackToBackExponential_sptr> peakfuncs, std::vector<bool> &vecfitgood,
132 std::vector<double> &vecchi2s);
136 std::vector<Functions::BackToBackExponential_sptr> peaks);
139 void setOverlappedPeaksConstraints(
const std::vector<Functions::BackToBackExponential_sptr> &peaks);
144 size_t maxiteration,
double &chi2);
161 std::vector<std::string> minimzernames, std::vector<size_t> maxiterations,
162 const std::vector<double> &dampfactors,
double &chi2);
167 const std::vector<Functions::BackToBackExponential_sptr> &peakfuncs,
168 const std::string &minimizername,
size_t maxiteration,
double &chi2);
171 void storeFunctionParameters(
const API::IFunction_sptr &function, std::map<std::string, double> ¶mmaps);
174 void restoreFunctionParameters(
const API::IFunction_sptr &function, std::map<std::string, double> parammap);
177 void calculatePeakFitBoundary(
size_t ileftpeak,
size_t irightpeak,
double &peakleftboundary,
178 double &peakrightboundary);
184 double ¢erleftbound,
double ¢errightbound,
int &errordirection);
190 std::pair<DataObjects::TableWorkspace_sptr, DataObjects::TableWorkspace_sptr> genPeakParametersWorkspace();
193 void cropWorkspace(
double tofmin,
double tofmax);
207 size_t m_wsIndex,
double &chi2);
212 double &peakleftbound,
double &peakrightbound);
227 double leftpeakbound,
double rightpeakbound);
238 double leftfwhm,
double rightfwhm,
double ¢er,
double &
sigma,
double &
height);
253 double getParameter(
const std::string &parname);
260 double leftbound,
double rightbound);
322 enum { ROBUSTFIT, TRUSTINPUTFIT } m_fitMode;
325 enum { HKLCALCULATION, FROMBRAGGTABLE } m_genPeakStartingValue;
349 double x = ((xf - x0) *
y - (xf * y0 - x0 * yf)) / (yf - y0);
356 double y = ((xf * y0 - x0 * yf) +
x * (yf - y0)) / (xf - x0);
363 size_t wsindexbkgd,
size_t wsindexpeak);
367 double &fwhm, std::string &errmsg);
Base class from which all concrete algorithm classes should be derived.
Implements FunctionDomain1D with its own storage in form of a std::vector.
FitPowderDiffPeaks : Fit peaks in powder diffraction pattern.
const std::string category() const override
Algorithm's category for identification overriding a virtual method.
double m_rightmostPeakLeftBound
Right most peak's left boundary.
double calculatePeakCentreTOF(int h, int k, int l)
Calculate a Bragg peak's centre in TOF from its Miller indices.
std::vector< double > m_chi2GoodFitPeaks
std::vector< int > m_rightmostPeakHKL
Right most peak HKL.
bool fitPeak(Functions::BackToBackExponential_sptr peak, Functions::BackgroundFunction_sptr background, double leftdev, double rightdev, size_t m_wsIndex, double &chi2)
Fit a single peak.
void calculate1PeakGroup(std::vector< size_t > peakindexes, Functions::BackgroundFunction_sptr background)
Calcualte the value of a single peak in a given range.
std::vector< int > m_minimumHKL
Minimum HKL.
double m_tofMin
TOF Min and TOF Max.
const std::vector< std::string > seeAlso() const override
Function to return all of the seeAlso (these are not validated) algorithms related to this algorithm....
double m_rightmostPeakRightBound
Right most peak's right boundary.
std::vector< bool > m_goodFit
Peak fitting status.
std::vector< double > m_peakFitChi2
Peak fitting information.
bool estimateFWHM(DataObjects::Workspace2D_sptr dataws, size_t wsindex, double tof_h, double &leftfwhm, double &rightfwhm)
Estimate FWHM for the peak observed.
bool fitSinglePeakConfidentX(Functions::BackToBackExponential_sptr peak)
Fit peak with confidence of the centre.
double m_minPeakHeight
Minimum peak height for peak to be refined.
std::pair< bool, double > doFitPeak_Old(DataObjects::Workspace2D_sptr dataws, Functions::BackToBackExponential_sptr peak, double guessedfwhm, bool calchi2)
Fit single peak without background.
bool m_fitPeakBackgroundComposite
Fit peak + background as the last step.
void fitPeaksGroup(std::vector< size_t > peakindexes)
Fit peaks in the same group (i.e., single peak or overlapped peaks)
int version() const override
Algorithm's version for identification overriding a virtual method.
DataObjects::TableWorkspace_sptr m_profileTable
Instrument profile parameter table workspace.
Geometry::UnitCell m_unitCell
Unit cell of powder crystal.
const std::string summary() const override
Summary of algorithms purpose.
std::vector< double > m_peakData
Data for each individual peaks. (HKL)^2, vector index, function values.
DataObjects::TableWorkspace_sptr m_peakParamTable
Bragg peak parameter.
void subtractBackground(DataObjects::Workspace2D_sptr dataws)
Estimate background.
bool findMaxHeight(API::MatrixWorkspace_sptr dataws, size_t wsindex, double xmin, double xmax, double ¢er, double ¢erleftbound, double ¢errightbound, int &errordirection)
Find max height (peak center)
std::vector< std::string > mPeakParameterNames
Peak parmeter names.
bool fitSinglePeakConfidentY(DataObjects::Workspace2D_sptr dataws, Functions::BackToBackExponential_sptr peak, double dampingfactor)
Fit peak with confident parameters.
bool estimateSinglePeakRange(Functions::BackToBackExponential_sptr peak, Functions::BackgroundFunction_sptr background, Functions::BackToBackExponential_sptr rightpeak, double fwhm, bool ismostright, size_t m_wsIndex, double &chi2)
Build a partial workspace from source.
std::map< std::string, double > m_instrumentParmaeters
Map for function (instrument parameter)
std::vector< std::pair< double, std::pair< std::vector< int >, Functions::BackToBackExponential_sptr > > > m_vecPeakFunctions
Sorted vector for peaks. double = d_h, vector = (HKL), peak.
int m_numPeaksLowerToMin
Number of peaks to fit lower to minimum HKL.
const std::string name() const override
Algorithm's name for identification overriding a virtual method.
bool m_confidentInInstrumentParameters
Flag to show whether input instrument parameters is trustful.
bool m_useGivenTOFh
Flag to use given Bragg peaks' centre in TOF.
int m_wsIndex
TOF vector of data workspace to process with.
API::MatrixWorkspace_sptr m_dataWS
Data.
std::vector< size_t > m_indexGoodFitPeaks
bool doFitBackground(DataObjects::Workspace2D_sptr dataws, Functions::BackgroundFunction_sptr background, double leftpeakbound, double rightpeakbound)
Fit background function by removing the peak properly.
Class to implement unit cell of crystals.
std::shared_ptr< IAlgorithm > IAlgorithm_sptr
shared pointer to Mantid::API::IAlgorithm
std::shared_ptr< ITableWorkspace > ITableWorkspace_sptr
shared pointer to Mantid::API::ITableWorkspace
std::shared_ptr< IFunction > IFunction_sptr
shared pointer to the function base class
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
std::shared_ptr< CompositeFunction > CompositeFunction_sptr
shared pointer to the composite function base class
double linearInterpolateX(double x0, double xf, double y0, double yf, double y)
Formular for linear iterpolation: X = [(xf-x0)*Y - (xf*y0-x0*yf)]/(yf-y0)
size_t findMaxValue(const std::vector< double > &Y)
Find maximum value.
bool observePeakParameters(const DataObjects::Workspace2D_sptr &dataws, size_t wsindex, double ¢re, double &height, double &fwhm, std::string &errmsg)
Estimate peak parameters;.
void estimateBackgroundCoarse(const DataObjects::Workspace2D_sptr &dataws, const Functions::BackgroundFunction_sptr &background, size_t wsindexraw, size_t wsindexbkgd, size_t wsindexpeak)
Estimate background for a pattern in a coarse mode.
double linearInterpolateY(double x0, double xf, double y0, double yf, double x)
Formula for linear interpolation: Y = ( (xf*y0-x0*yf) + x*(yf-y0) )/(xf-x0)
std::string getFunctionInfo(const API::IFunction_sptr &function)
Get function parameter name, value and etc information in string.
std::shared_ptr< BackToBackExponential > BackToBackExponential_sptr
std::shared_ptr< BackgroundFunction > BackgroundFunction_sptr
std::shared_ptr< Workspace2D > Workspace2D_sptr
shared pointer to Mantid::DataObjects::Workspace2D
std::shared_ptr< TableWorkspace > TableWorkspace_sptr
shared pointer to Mantid::DataObjects::TableWorkspace
Helper class which provides the Collimation Length for SANS instruments.