Mantid
Loading...
Searching...
No Matches
ThermalNeutronDtoTOFFunction.h
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2012 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
7#pragma once
8
12#include "MantidAPI/Jacobian.h"
14#include "MantidCurveFitting/DllConfig.h"
15#include "MantidKernel/Logger.h"
16#include "MantidKernel/System.h"
17
18#include <cmath>
19#include <gsl/gsl_sf_erf.h>
20
21namespace Mantid {
22namespace CurveFitting {
23namespace Functions {
24
27class MANTID_CURVEFITTING_DLL ThermalNeutronDtoTOFFunction : virtual public API::IFunction1D,
28 public API::ParamFunction {
29public:
31 void function1D(double *out, const double *xValues, const size_t nData) const override;
32
34 std::string name() const override { return "ThermalNeutronDtoTOFFunction"; }
35
37 const std::string category() const override { return "General"; }
38
40 void function1D(std::vector<double> &out, const std::vector<double> &xValues) const;
41
42protected:
44 void init() override;
45
46private:
48 inline double corefunction(double dh, double dtt1, double dtt1t, double dtt2t, double zero, double zerot,
49 double width, double tcross) const;
50
52 void functionDerivLocal(API::Jacobian *, const double *, const size_t);
53
55 // void functionDeriv(const API::FunctionDomain& domain, API::Jacobian&
56 // jacobian);
57
59 void functionDeriv1D(API::Jacobian *out, const double *xValues, const size_t nData) override;
60
61 // void functionDeriv(const API::FunctionDomain& domain, API::Jacobian&
62 // jacobian);
63};
64
65using ThermalNeutronDtoTOFFunction_sptr = std::shared_ptr<ThermalNeutronDtoTOFFunction>;
66
68inline double calThermalNeutronTOF(double dh, double dtt1, double dtt1t, double dtt2t, double zero, double zerot,
69 double width, double tcross) {
70 double n = 0.5 * gsl_sf_erfc(width * (tcross - 1 / dh));
71 double Th_e = zero + dtt1 * dh;
72 double Th_t = zerot + dtt1t * dh - dtt2t / dh;
73 double tof_h = n * Th_e + (1 - n) * Th_t;
74
75 return tof_h;
76}
77
78} // namespace Functions
79} // namespace CurveFitting
80} // namespace Mantid
This is a specialization of IFunction for functions of one real argument.
Definition: IFunction1D.h:43
Represents the Jacobian in IFitFunction::functionDeriv.
Definition: Jacobian.h:22
Implements the part of IFunction interface dealing with parameters.
Definition: ParamFunction.h:33
const std::string category() const override
Overwrite IFunction.
std::string name() const override
overwrite IFunction base class methods
double corefunction(double dh, double dtt1, double dtt1t, double dtt2t, double zero, double zerot, double width, double tcross) const
Core function (inline) to calcualte TOF_h from d-spacing.
double calThermalNeutronTOF(double dh, double dtt1, double dtt1t, double dtt2t, double zero, double zerot, double width, double tcross)
Calcualte TOF from d-spacing value for thermal neutron.
std::shared_ptr< ThermalNeutronDtoTOFFunction > ThermalNeutronDtoTOFFunction_sptr
Helper class which provides the Collimation Length for SANS instruments.