LeBailFit.h

Go to the documentation of this file.

// Mantid Repository : https://github.com/mantidproject/mantid
//
// Copyright © 2012 ISIS Rutherford Appleton Laboratory UKRI,
//   NScD Oak Ridge National Laboratory, European Spallation Source,
//   Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
// SPDX - License - Identifier: GPL - 3.0 +
#pragma once
 
#include "MantidAPI/Algorithm.h"
#include "MantidAPI/CompositeFunction.h"
#include "MantidAPI/IFunction.h"
#include "MantidAPI/ITableWorkspace_fwd.h"
#include "MantidAPI/MatrixWorkspace_fwd.h"
#include "MantidCurveFitting/Algorithms/LeBailFunction.h"
#include "MantidCurveFitting/DllConfig.h"
#include "MantidCurveFitting/Functions/BackgroundFunction.h"
#include "MantidCurveFitting/Functions/ThermalNeutronBk2BkExpConvPVoigt.h"
#include "MantidDataObjects/TableWorkspace.h"
#include "MantidDataObjects/Workspace2D.h"
 
#include <gsl/gsl_sf_erf.h>
 
#include <span>
 
namespace Mantid {
namespace HistogramData {
class HistogramX;
class HistogramY;
} // namespace HistogramData
 
namespace CurveFitting {
namespace Algorithms {
 
struct Parameter {
  // Regular
  std::string name;
  double curvalue = 0;
  double prevalue = 0;
  double minvalue = 0;
  double maxvalue = 0;
  bool fit = false;
  double stepsize = 0;
  double fiterror = 0;
  // Monte Carlo
  bool nonnegative = false;
  double mcA0 = 0;
  double mcA1 = 0;
  // Monte Carlo record
  double sumstepsize = 0;
  double maxabsstepsize = 0;
  double maxrecordvalue = 0;
  double minrecordvalue = 0;
  size_t numpositivemove = 0;
  size_t numnegativemove = 0;
  size_t numnomove = 0;
  int movedirection = 0;
};
 
class MANTID_CURVEFITTING_DLL LeBailFit final : public API::Algorithm {
public:
  enum FunctionMode { CALCULATION, FIT, BACKGROUNDPROCESS, MONTECARLO };
 
  LeBailFit();
 
  const std::string name() const override { return "LeBailFit"; }
  const std::string summary() const override { return "Do LeBail Fit to a spectrum of powder diffraction data. "; }
 
  int version() const override { return 1; }
  const std::vector<std::string> seeAlso() const override { return {"CreateLeBailFitInput", "FitPowderDiffPeaks"}; }
 
  const std::string category() const override { return "Diffraction\\Fitting"; }
 
private:
  // Implement abstract Algorithm methods
  void init() override;
  // Implement abstract Algorithm methods
  void exec() override;
 
  void processInputProperties();
 
  //--------------  Pattern Calculation & Minimizing  -------------------
  void execPatternCalculation();
 
  void execRefineBackground();
 
  //--------------  Functions to set up the Le Bail Fit -----------------
  void createLeBailFunction();
 
  API::MatrixWorkspace_sptr cropWorkspace(const API::MatrixWorkspace_sptr &inpws, size_t wsindex);
 
  void processInputBackground();
 
  //--------------  Le Bail Formular: Calculate Peak Intensities ------------
  bool calculatePeaksIntensities(API::MatrixWorkspace_sptr dataws, size_t workspaceindex, bool zerobackground,
                                 std::vector<double> &allpeaksvalues);
 
  //--------------  Import and Export ---------------------------------------
  void parseInstrumentParametersTable();
 
  void parseBraggPeaksParametersTable();
 
  void parseBackgroundTableWorkspace(const DataObjects::TableWorkspace_sptr &bkgdparamws,
                                     std::vector<std::string> &bkgdparnames, std::vector<double> &bkgdorderparams);
 
  void exportBraggPeakParameterToTable();
 
  void exportInstrumentParameterToTable(const std::map<std::string, Parameter> &parammap);
 
  void createOutputDataWorkspace();
 
  //--------------  Random Walk Suite ----------------------------------------
  void execRandomWalkMinimizer(size_t maxcycles, std::map<std::string, Parameter> &parammap);
 
  void doMarkovChain(const std::map<std::string, Parameter> &parammap, const Mantid::HistogramData::HistogramX &vecX,
                     const Mantid::HistogramData::HistogramY &vecPurePeak,
                     const Mantid::HistogramData::HistogramY &vecBkgd, size_t maxcycles, const Kernel::Rfactor &startR,
                     int randomseed);
 
  void setupBuiltInRandomWalkStrategy();
 
  void setupRandomWalkStrategyFromTable(const DataObjects::TableWorkspace_sptr &tablews);
 
  void addParameterToMCMinimize(std::vector<std::string> &parnamesforMC, const std::string &parname);
 
  bool calculateDiffractionPattern(const Mantid::HistogramData::HistogramX &vecX,
                                   const Mantid::HistogramData::HistogramY &vecY, bool inputraw, bool outputwithbkgd,
                                   const Mantid::HistogramData::HistogramY &vecBkgd, std::vector<double> &values,
                                   Kernel::Rfactor &rfactor);
 
  bool acceptOrDeny(Kernel::Rfactor currR, Kernel::Rfactor newR);
 
  bool proposeNewValues(const std::vector<std::string> &mcgroup, Kernel::Rfactor r,
                        std::map<std::string, Parameter> &curparammap, std::map<std::string, Parameter> &newparammap,
                        bool prevBetterRwp);
 
  double limitProposedValueInBound(const Parameter &param, double newvalue, double direction, int choice);
 
  void bookKeepBestMCResult(const std::map<std::string, Parameter> &parammap, std::span<double const> bkgddata,
                            Kernel::Rfactor rfactor, size_t istep);
 
  void applyParameterValues(const std::map<std::string, Parameter> &srcparammap,
                            std::map<std::string, Parameter> &tgtparammap);
 
  void storeBackgroundParameters(std::vector<double> &bkgdparamvec);
 
  void recoverBackgroundParameters(const std::vector<double> &bkgdparamvec);
 
  void proposeNewBackgroundValues();
 
  //--------------------------------------------------------------------------------------------
 
  LeBailFunction_sptr m_lebailFunction;
 
  API::MatrixWorkspace_sptr m_dataWS;
  DataObjects::Workspace2D_sptr m_outputWS;
  DataObjects::TableWorkspace_sptr parameterWS;
  DataObjects::TableWorkspace_sptr reflectionWS;
 
  size_t m_wsIndex;
 
  double m_startX, m_endX;
 
  std::vector<std::pair<std::vector<int>, double>> m_inputPeakInfoVec;
 
  Functions::BackgroundFunction_sptr m_backgroundFunction;
 
  std::map<std::string, Parameter> m_funcParameters; // char = f: fit... = t: tie to value
  std::map<std::string, double> m_origFuncParameters;
 
  std::map<std::string, double> convertToDoubleMap(std::map<std::string, Parameter> &inmap);
 
 
  std::string m_peakType;
 
  std::string m_backgroundType;
 
  std::vector<double> m_backgroundParameters;
  std::vector<std::string> m_backgroundParameterNames;
  unsigned int m_bkgdorder;
 
  int mPeakRadius;
 
  double m_lebailFitChi2;
  double m_lebailCalChi2;
 
  std::string mMinimizer;
  double m_dampingFactor;
 
  bool m_inputParameterPhysical;
 
  FunctionMode m_fitMode;
 
  double m_indicatePeakHeight;
 
  //-------------------------- Monte Carlo Variables--------------------------
  std::map<int, std::vector<std::string>> m_MCGroups;
  size_t m_numMCGroups;
 
  double m_bestRwp;
  double m_bestRp;
 
  std::map<std::string, Parameter> m_bestParameters;
  std::vector<double> m_bestBackgroundData;
  size_t m_bestMCStep;
 
  size_t m_numMinimizeSteps;
 
  double m_Temperature;
 
  bool m_useAnnealing;
 
  enum { RANDOMWALK, DRUNKENWALK } m_walkStyle;
 
  double m_minimumPeakHeight;
 
  bool m_tolerateInputDupHKL2Peaks;
 
  //------------------------ Background Refinement Variables
  //-----------------------
  std::vector<std::string> m_bkgdParameterNames;
  size_t m_numberBkgdParameters;
  std::vector<double> m_bkgdParameterBuffer;
  std::vector<double> m_bestBkgdParams;
  int m_roundBkgd;
  std::vector<double> m_bkgdParameterStepVec;
 
  double m_peakCentreTol;
};
 
void writeRfactorsToFile(const std::vector<double> &vecX, const std::vector<Kernel::Rfactor> &vecR,
                         const std::string &filename);
 
} // namespace Algorithms
} // namespace CurveFitting
} // namespace Mantid