diff --git a/configure.ac b/configure.ac index 933aa13e..9efd1188 100644 --- a/configure.ac +++ b/configure.ac @@ -1288,6 +1288,7 @@ AC_CONFIG_FILES([Makefile \ src/external/libPhotoMeissner/classes/Makefile \ src/external/libGbGLF/Makefile \ src/external/libBNMR/Makefile \ + src/external/libBNMR/libLineProfile/Makefile \ src/musredit_qt5/Makefile \ src/musredit/Makefile \ src/musrgui/Makefile \ diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index e8ba8f90..23ce8455 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -60,6 +60,7 @@ add_executable(dump_header git-revision.h dump_header.cpp) target_compile_options(dump_header BEFORE PRIVATE "-DHAVE_CONFIG_H") target_include_directories(dump_header BEFORE PRIVATE + $ $ $ $ @@ -77,6 +78,7 @@ add_executable(msr2data git-revision.h msr2data.cpp) target_compile_options(msr2data BEFORE PRIVATE "-DHAVE_CONFIG_H") target_include_directories(msr2data BEFORE PRIVATE + $ $ $ $ diff --git a/src/external/libBNMR/CMakeLists.txt b/src/external/libBNMR/CMakeLists.txt index 1f10f260..356b61f5 100644 --- a/src/external/libBNMR/CMakeLists.txt +++ b/src/external/libBNMR/CMakeLists.txt @@ -1,5 +1,8 @@ # - beta-NMR library ---------------------------------------------------------- +#--- add subdirectory libLineProfile +add_subdirectory(libLineProfile) + #--- generate necessary dictionaries ------------------------------------------ set(MUSRFIT_INC ${CMAKE_SOURCE_DIR}/src/include) diff --git a/src/external/libBNMR/Makefile.am b/src/external/libBNMR/Makefile.am index f099256c..750017fe 100644 --- a/src/external/libBNMR/Makefile.am +++ b/src/external/libBNMR/Makefile.am @@ -1,4 +1,5 @@ ## Process this file with automake to create Makefile.in +SUBDIRS = libLineProfile h_sources = \ TBNMR.h diff --git a/src/external/libBNMR/TBNMRLinkDef.h b/src/external/libBNMR/TBNMRLinkDef.h index ab23c35f..3eff15c5 100644 --- a/src/external/libBNMR/TBNMRLinkDef.h +++ b/src/external/libBNMR/TBNMRLinkDef.h @@ -35,9 +35,7 @@ #pragma link off all classes; #pragma link off all functions; -#pragma link C++ class TBNMR+; #pragma link C++ class ExpRlx+; #pragma link C++ class SExpRlx+; -#pragma link C++ class MLRes+; #endif //__CINT__ diff --git a/src/external/libBNMR/libLineProfile/CMakeLists.txt b/src/external/libBNMR/libLineProfile/CMakeLists.txt new file mode 100644 index 00000000..3663db0e --- /dev/null +++ b/src/external/libBNMR/libLineProfile/CMakeLists.txt @@ -0,0 +1,55 @@ +# - beta-NMR LineProfile library ---------------------------------------------------------- + +#--- generate necessary dictionaries ------------------------------------------ +set(MUSRFIT_INC ${CMAKE_SOURCE_DIR}/src/include) + +root_generate_dictionary( + libLineProfileDict + -I${FFTW3_INCLUDE_DIR} + -I${MUSRFIT_INC} + libLineProfile.h + LINKDEF libLineProfileLinkDef.h + MODULE libLineProfile +) + +#--- lib creation ------------------------------------------------------------- +add_library(LineProfile SHARED + libLineProfile.cpp + libLineProfileDict.cxx +) + +#--- make sure that the include directory is found ---------------------------- +target_include_directories( + LineProfile BEFORE PRIVATE + $ + $ + $ +) + +#--- set target properties, e.g. version -------------------------------------- +set_target_properties(LineProfile + PROPERTIES + VERSION "1.0.0" +) + +#--- add library dependencies ------------------------------------------------- +target_link_libraries(LineProfile ${ROOT_LIBRARIES} PUserFcnBase) + +#--- install libLineProfile solib ---------------------------------------------------- +install(TARGETS LineProfile DESTINATION lib) + +#--- install root pcm's and rootmaps ------------------------------------------ +install( + FILES ${CMAKE_CURRENT_BINARY_DIR}/libLineProfile_rdict.pcm + ${CMAKE_CURRENT_BINARY_DIR}/libLineProfile.rootmap + DESTINATION lib +) + +#--- install libLineProfile header --------------------------------------------------- +install( + FILES + libLineProfile.h + DESTINATION + include +) + diff --git a/src/external/libBNMR/libLineProfile/LineProfile.pdf b/src/external/libBNMR/libLineProfile/LineProfile.pdf new file mode 100644 index 00000000..5827e620 Binary files /dev/null and b/src/external/libBNMR/libLineProfile/LineProfile.pdf differ diff --git a/src/external/libBNMR/libLineProfile/LineProfile.tex b/src/external/libBNMR/libLineProfile/LineProfile.tex new file mode 100644 index 00000000..f655496d --- /dev/null +++ b/src/external/libBNMR/libLineProfile/LineProfile.tex @@ -0,0 +1,190 @@ +\documentclass[twoside]{article} + +\usepackage[english]{babel} +%\usepackage{a4} +\usepackage{amssymb,amsmath,bm} +\usepackage{graphicx,tabularx} +\usepackage{fancyhdr} +\usepackage{array} +\usepackage{float} +\usepackage{hyperref} +\usepackage{xspace} +\usepackage{rotating} +\usepackage{dcolumn} +\usepackage{geometry} +\usepackage{color} + +\geometry{a4paper,left=20mm,right=20mm,top=20mm,bottom=20mm} + +% \setlength{\topmargin}{10mm} +% \setlength{\topmargin}{-13mm} +% % \setlength{\oddsidemargin}{0.5cm} +% % \setlength{\evensidemargin}{0cm} +% \setlength{\oddsidemargin}{1cm} +% \setlength{\evensidemargin}{1cm} +% \setlength{\textwidth}{15cm} +\setlength{\textheight}{23.8cm} + +\pagestyle{fancyplain} +\addtolength{\headwidth}{0.6cm} +\fancyhead{}% +\fancyhead[RE,LO]{\bf \textsc{LineProfile}}% +\fancyhead[LE,RO]{\thepage} +\cfoot{--- J.~A.~Krieger -- \today~ ---} +\rfoot{\includegraphics[width=2cm]{PSI-Logo_narrow.jpg}} + +\DeclareMathAlphabet{\bi}{OML}{cmm}{b}{it} + +\newcommand{\mean}[1]{\langle #1 \rangle} +\newcommand{\lem}{LE-$\mu$SR\xspace} +\newcommand{\lemhead}{LE-$\bm{\mu}$SR\xspace} +\newcommand{\musr}{$\mu$SR\xspace} +\newcommand{\musrhead}{$\bm{\mu}$SR\xspace} +\newcommand{\trimsp}{\textsc{trim.sp}\xspace} +\newcommand{\musrfithead}{MUSRFIT\xspace} +\newcommand{\musrfit}{\textsc{musrfit}\xspace} +\newcommand{\gapint}{\textsc{GapIntegrals}\xspace} +\newcommand{\YBCO}{YBa$_{2}$Cu$_{3}$O$_{7-\delta}$\xspace} +\newcommand{\YBCOhead}{YBa$_{\bm{2}}$Cu$_{\bm{3}}$O$_{\bm{7-\delta}}$\xspace} + +\newcolumntype{d}[1]{D{.}{.}{#1}} +\newcolumntype{C}[1]{>{\centering\arraybackslash}p{#1}} + +\begin{document} +% Header info -------------------------------------------------- +\thispagestyle{empty} +\noindent +\begin{tabular}{@{\hspace{-0.2cm}}l@{\hspace{6cm}}r} +\noindent\includegraphics[width=3.4cm]{PSI-Logo_narrow.jpg} & + {\Huge\sf Memorandum} +\end{tabular} +% +\vskip 1cm +% +\begin{tabular}{@{\hspace{-0.5cm}}ll@{\hspace{4cm}}ll} +Date: & \today & & \\[3ex] +From: & J.~A.~Krieger & \\ +E-Mail: & \verb?jonas.krieger@psi.ch? && +\end{tabular} +% +\vskip 0.3cm +\noindent\hrulefill +\vskip 1cm +% +\section*{\musrfithead plug-in for simple $\beta$-NMR resonance line shapes}% +This library contains useful functions to fit NMR and $\beta$-NMR line shapes. +The functional form of the powder averages was taken from +\href{http://dx.doi.org/10.1007/978-3-642-68756-3_2}{M. Mehring, Principles +of High Resolution NMR in Solids (Springer 1983)}. +% +The \texttt{libLineProfile} library currently contains the following functions: +\begin{description} + \item[LineGauss] + \begin{equation} + A(f)=e^{-\frac{4\ln 2 (f-f_0)^2}{ \sigma^2}} + \end{equation} + Gaussian line shape around $f_0$ with width $\sigma$ and height~$1$.\\[1.5ex] + \musrfit theory line: \verb?userFcn libLineProfile LineGauss 1 2?\\[1.5ex] + Parameters: $f_0$, $\sigma$. + \item[LineLaplace] + \begin{equation} + A(f)=e^{-2\ln 2 \left|\frac{f-f_0}{\sigma}\right|} + \end{equation} + Laplaceian line shape around $f_0$ with width $\sigma$ and +height~$1$.\\[1.5ex] + \musrfit theory line: \verb?userFcn libLineProfile LineLaplace 1 2? +\\[1.5ex] + Parameters: $f_0$, $\sigma$. + \item[LineLorentzian] + \begin{equation} + A(f)= +\frac{\sigma^2}{4(f-f_0)^2+\sigma^2} + \end{equation} + Lorentzian line shape around $f_0$ with width $\sigma$ and +height~$1$.\\[1.5ex] + \musrfit theory line: \verb?userFcn libLineProfile LineLorentzian 1 2? +\\[1.5ex] + Parameters: $f_0$, $\sigma$. + \item[LineSkewLorentzian] + \begin{equation} + A(f)= \frac{\sigma*\sigma_a}{4(f-f_0)^2+\sigma_a^2}, \quad \sigma_a=\frac{2\sigma}{1+e^{a(f-f_0)}} + \end{equation} + Skewed Lorentzian line shape around $f_0$ with width $\sigma$, +height~$1$ and skewness parameter $a$.\\[1.5ex] + \musrfit theory line: \verb?userFcn libLineProfile LineSkewLorentzian 1 2 3? +\\[1.5ex] + Parameters: $f_0$, $\sigma$, $a$. + \item[LineSkewLorentzian2] + \begin{equation} + A(f)= \left\{\begin{matrix}\frac{{\sigma_1}^2}{4{(f-f_0)}^2+{\sigma_1}^2},&ff_0\end{matrix}\right. + \end{equation} + Skewed Lorentzian line shape around $f_0$ with height~$1$ and widths $\sigma_1$, + and $\sigma_2$.\\[1.5ex] + \musrfit theory line: \verb?userFcn libLineProfile LineSkewLorentzian2 1 2 3? +\\[1.5ex] + Parameters: $f_0$, $\sigma_1$, $\sigma_2$. + + +\item[PowderLineAxialLor] + \begin{equation} + A(f)= I_{\mathrm ax}(f)\circledast\left( \frac{\sigma^2}{4f^2+\sigma^2} \right) + \end{equation} + Powder average of a axially symmetric interaction, convoluted with a Lorentzian. + \begin{equation}\label{eq:Iax} + I_{\mathrm ax}(f)=\left\{\begin{matrix} \frac{1}{2\sqrt{(f_\parallel-f_\perp)(f-f_\perp)}}& f\in(f_\perp,f_\parallel)\cup(f_\parallel,f_\perp)\\[6pt] 0 & \text{otherwise}\end{matrix} \right. + \end{equation} + The maximal height of the curve is normalized to $\sim$1. + \\[1.5ex] + \musrfit theory line: \verb?userFcn libLineProfile PowderLineAxialLor 1 2 3? +\\[1.5ex] + Parameters: $f_\parallel$, $f_\perp$, $\sigma$. + +\item[PowderLineAxialGss] + \begin{equation} + A(f)= I_{\mathrm ax}(f)\circledast\left(e^{-\frac{4\ln 2 (f-f_0)^2}{ +\sigma^2}} \right) + \end{equation} + Powder average of a axially symmetric interaction (Eq.~\ref{eq:Iax}), convoluted with a Gaussian. The maximal height of the curve is normalized to $\sim$1. + \\[1.5ex] + \musrfit theory line: \verb?userFcn libLineProfile PowderLineAxialGss 1 2 3? +\\[1.5ex] + Parameters: $f_\parallel$, $f_\perp$, $\sigma$. + +\item[PowderLineAsymLor] + \begin{equation} + A(f)= I(f)\circledast\left( \frac{\sigma^2}{4f^2+\sigma^2} \right) + \end{equation} + Powder average of a asymmetric interaction, convoluted with a Lorentzian. + Assume without loss of generality that $f_1f_2 \\[9pt] + \frac{K(m)}{\pi\sqrt{(f_3-f)(f_2-f_1)}},& f_2>f\geq f_1\\[9pt] + 0 & \text{otherwise} + \end{matrix} \right. \\ + m&=\left\{\begin{matrix} + \frac{(f_2-f_1)(f_3-f)}{(f_3-f_2)(f-f_1)},& f_3\geq f>f_2 \\[6pt] + \frac{(f-f_1)(f_3-f_2)}{(f_3-f)(f_2-f_1)},& f_2>f\geq f_1\\[6pt] + \end{matrix} \right. \\\label{eq:Kofm} + K(m)&=\int_0^{\pi/2}\frac{\mathrm d\varphi}{\sqrt{1-m^2\sin^2{\varphi}}}, + \end{align} + where $K(m)$ is the complete elliptic integral of the first kind. + Note that $f_1 Building shared library $(SHLIB) ..." + /bin/rm -f $(SHLIB) + $(LD) $(OBJS) $(SOFLAGS) -o $(SHLIB) + @echo "done" + +# clean up: remove all object file (and core files) +# semicolon needed to tell make there is no source +# for this target! +# +clean:; @rm -f $(OBJS) *Dict* core* + @echo "---> removing $(OBJS)" + +# +$(OBJS): %.o: %.cpp + $(CXX) $(INCLUDES) $(CXXFLAGS) -c $< + +# Generate the ROOT CINT dictionary + +LineProfileDict.cpp: LineProfile.h LineProfileLinkDef.h + @echo "Generating dictionary $@..." + rootcint -f $@ -c -p -I$(ROOTINCLUDE) $^ + +install: all + @echo "Installing shared lib: libLineProfile.so" +ifeq ($(OS),LINUX) + cp -pv $(SHLIB) $(ROOTSYS)/lib +endif + +uninstall:; +ifeq ($(OS),LINUX) + rm $(ROOTSYS)/lib/$(SHLIB) +endif + @echo "Installing shared lib: libLineProfile.so" + diff --git a/src/external/libBNMR/libLineProfile/PSI-Logo_narrow.jpg b/src/external/libBNMR/libLineProfile/PSI-Logo_narrow.jpg new file mode 100644 index 00000000..9bda1610 Binary files /dev/null and b/src/external/libBNMR/libLineProfile/PSI-Logo_narrow.jpg differ diff --git a/src/external/libBNMR/libLineProfile/libLineProfile.cpp b/src/external/libBNMR/libLineProfile/libLineProfile.cpp new file mode 100644 index 00000000..48269abb --- /dev/null +++ b/src/external/libBNMR/libLineProfile/libLineProfile.cpp @@ -0,0 +1,847 @@ +/**************************************************************************/ + +/*************************************************************************** + * Copyright (C) 2017 by Jonas A. Krieger * + * jonas.krieger@psi.ch * + * * + * This program is free software; you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation; either version 2 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program; if not, write to the * + * Free Software Foundation, Inc., * + * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. * +***************************************************************************/ + +#include "libLineProfile.h" +#include //for testing purposes + +//Implement helperfunctions + +Double_t GaussianShape(Double_t x, Double_t position, Double_t width) { + if(!width){ // width=0 + if(x==position) return 1.0; + return 0.0;} + return exp(-2.7725887222397811*pow((x-position)/width,2)); +} + +Double_t LaplacianShape(Double_t x, Double_t position, Double_t width) { + if(!width){ // width=0 + if(x==position) return 1.0; + return 0.0;} + return exp(-1.3862943611198906*abs((x-position)/width)); +} + +Double_t LorentzianShape(Double_t x, Double_t position, Double_t width) { + if(!width){ // width=0 + if(x==position) return 1.0; + return 0.0;} + + // return 2/PI*(width/(4*pow(x-position,2)+pow(width,2)));//constant Area + return (1/(4*(x-position)*(x-position)/width/width+1));//constant Height +} + + +Double_t IAxial(Double_t x, Double_t omega_par, Double_t omega_per){ + if(omega_par==omega_per) return 0; // delta function at x==omega_per, taken care of in the userfunction + if((x<=min(omega_par,omega_per))||(x>=max(omega_par,omega_per))) return 0; + if (abs(x-omega_per)<0.001*abs(omega_par-omega_per)){return 1./2/pow((omega_par-omega_per)*(0.03*(omega_par-omega_per)),0.5);} + return 1./2/pow((omega_par-omega_per)*(x-omega_per),0.5); +} + + +Double_t IAsym(Double_t x, Double_t omega_center, Double_t omega_min,Double_t omega_max){ + if (omega_center==omega_min&&omega_min==omega_max) return 0; // delta function at x==omega_per, taken care of in the userfunction + Double_t cutoff=0.001; + if (abs(omega_center-omega_min)=omega_max)) return 0; + if (abs(omega_center-x)0){//omega_center &par) const { + assert(par.size()==2); // make sure the number of parameters handed to the function is correct + Double_t position=par[0]; + Double_t width=par[1]; + + return GaussianShape( x, position, width); +} + +ClassImp(LineLaplace) // for the ROOT dictionary + +Double_t LineLaplace::operator()(Double_t x, const vector &par) const { + assert(par.size()==2); // make sure the number of parameters handed to the function is correct + Double_t position=par[0]; + Double_t width=par[1]; + return LaplacianShape( x, position, width); + +} + + +ClassImp(LineLorentzian) // for the ROOT dictionary + +Double_t LineLorentzian::operator()(Double_t x, const vector &par) const { + assert(par.size()==2); // make sure the number of parameters handed to the function is correct + Double_t position=par[0]; + Double_t width=par[1]; + return LorentzianShape( x, position, width); + +} + +ClassImp(LineSkewLorentzian) // for the ROOT dictionary + +Double_t LineSkewLorentzian::operator()(Double_t x, const vector &par) const { + assert(par.size()==3); // make sure the number of parameters handed to the function is correct + Double_t position=par[0]; + Double_t width=par[1]; + Double_t a=par[2]; + + + // sigmoidal variation of width, Stancik2008 dx.doi.org/10.1016/j.vibspec.2008.02.009 + Double_t skewwidth=2*width/(1+exp(a*(x-position))); + + return width/skewwidth*LorentzianShape( x, position, skewwidth); +} + + +ClassImp(LineSkewLorentzian2) // for the ROOT dictionary + +Double_t LineSkewLorentzian2::operator()(Double_t x, const vector &par) const { + assert(par.size()==3); // make sure the number of parameters handed to the function is correct + Double_t position=par[0]; + Double_t widthlow=par[1]; //width below position + Double_t widthhigh=par[2]; //width above position + + if(x==position){return 1.0;} + else if(x ¶m){ + +assert(param.size() == 3); + +// check that param are new and hence a calculation is needed (except the phase, which is not needed here) +Bool_t newParams = false; +if (fPrevParam.size() == 0) { + for (UInt_t i=0; i0 +if(width<0) {width=-width;} + +Double_t omin=min(omega_par,omega_per); +Double_t omax=max(omega_par,omega_per); +Double_t x=omega_per; + +//Lorentzian window: maximally +- 15*width +//Lineshape Window [omin omax] +//Determine necessary convolution window +Double_t convmin=-15*width; +Double_t convmax=-convmin; +if(x-convmax>omin){convmax=x-omin;} +if(x-convmin0)) fNormalization=1; +return; +} + + +/******************************************************/ +// Implementation PowderLineAxialLor +ClassImp(PowderLineAxialLor) // for the ROOT dictionary + +PowderLineAxialLor::PowderLineAxialLor(){ + fValid = false; + fInvokedGlobal = false; + fGlobalUserFcn = 0; +} + + +PowderLineAxialLor::~PowderLineAxialLor() +{ + if ((fGlobalUserFcn != 0) && fInvokedGlobal) { + delete fGlobalUserFcn; + fGlobalUserFcn = 0; + } +} + + + +void PowderLineAxialLor::SetGlobalPart(vector &globalPart, UInt_t idx) +{ + fIdxGlobal = static_cast(idx); + + if ((Int_t)globalPart.size() <= fIdxGlobal) { // global user function not present, invoke it + fGlobalUserFcn = new PowderLineAxialLorGlobal(); + if (fGlobalUserFcn == 0) { // global user function object couldn't be invoked -> error + fValid = false; + cerr << endl << ">> PowderLineAxialLor::SetGlobalPart(): **ERROR** Couldn't invoke global user function object, sorry ..." << endl; + } else { // global user function object could be invoked -> resize to global user function vector and keep the pointer to the corresponding object + globalPart.resize(fIdxGlobal+1); + globalPart[fIdxGlobal] = dynamic_cast(fGlobalUserFcn); + fValid = true; + fInvokedGlobal = true; + } + } else { // global user function already present hence just retrieve a pointer to it + fValid = true; + fGlobalUserFcn = (PowderLineAxialLorGlobal*)globalPart[fIdxGlobal]; + } +} + + + +Double_t PowderLineAxialLor::operator()(Double_t x, const vector &par) const { + assert(par.size()==3); // make sure the number of parameters handed to the function is correct + Double_t omega_par=par[0]; + Double_t omega_per=par[1]; //middle of the line + Double_t width=par[2];//width of the Lorentzian + + + //ensure width>0 + if(width<0) width=-width; + + //calculate normalization globally + fGlobalUserFcn->CalcNormalization(par); + + // extract Normalization from the global object + Double_t Norm = fGlobalUserFcn->GetNormalization(); + + //symmetric situation has a delta function: + if(abs(omega_par-omega_per)<0.01*abs(width)){ + cout<omin){convmax=x-omin;} + if(x-convminGetSteps(); + for(Double_t y=convmin;y ¶m){ + +assert(param.size() == 3); + +// check that param are new and hence a calculation is needed (except the phase, which is not needed here) +Bool_t newParams = false; +if (fPrevParam.size() == 0) { + for (UInt_t i=0; i0 +if(width<0) width=-width; + + +Double_t omin=min(omega_par,omega_per); +Double_t omax=max(omega_par,omega_per); +Double_t x=omega_per; + +//Gaussian window: maximally +- 15*width +//Lineshape Window [omin omax] +//Determine necessary convolution window +Double_t convmin=-15*width; +Double_t convmax=-convmin; +if(x-convmax>omin){convmax=x-omin;} +if(x-convmin0)) fNormalization=1; +return; +} + + +/******************************************************/ +// Implementation PowderLineAxialGss +ClassImp(PowderLineAxialGss) // for the ROOT dictionary + +PowderLineAxialGss::PowderLineAxialGss(){ + fValid = false; + fInvokedGlobal = false; + fGlobalUserFcn = 0; +} + + +PowderLineAxialGss::~PowderLineAxialGss() +{ + if ((fGlobalUserFcn != 0) && fInvokedGlobal) { + delete fGlobalUserFcn; + fGlobalUserFcn = 0; + } +} + + + +void PowderLineAxialGss::SetGlobalPart(vector &globalPart, UInt_t idx) +{ + fIdxGlobal = static_cast(idx); + + if ((Int_t)globalPart.size() <= fIdxGlobal) { // global user function not present, invoke it + fGlobalUserFcn = new PowderLineAxialGssGlobal(); + if (fGlobalUserFcn == 0) { // global user function object couldn't be invoked -> error + fValid = false; + cerr << endl << ">> PowderLineAxialGss::SetGlobalPart(): **ERROR** Couldn't invoke global user function object, sorry ..." << endl; + } else { // global user function object could be invoked -> resize to global user function vector and keep the pointer to the corresponding object + globalPart.resize(fIdxGlobal+1); + globalPart[fIdxGlobal] = dynamic_cast(fGlobalUserFcn); + fValid = true; + fInvokedGlobal = true; + } + } else { // global user function already present hence just retrieve a pointer to it + fValid = true; + fGlobalUserFcn = (PowderLineAxialGssGlobal*)globalPart[fIdxGlobal]; + } +} + + + +Double_t PowderLineAxialGss::operator()(Double_t x, const vector &par) const { + assert(par.size()==3); // make sure the number of parameters handed to the function is correct + Double_t omega_par=par[0]; + Double_t omega_per=par[1]; //middle of the line + Double_t width=par[2];//width of the Gaussian + + + //ensure width>0 + if(width<0) width=-width; + + //calculate normalization globally + fGlobalUserFcn->CalcNormalization(par); + + // extract Normalization from the global object + Double_t Norm = fGlobalUserFcn->GetNormalization(); + + //symmetric situation has a delta function: + if(abs(omega_par-omega_per)<0.01*abs(width)){ + cout<omin){convmax=x-omin;} + if(x-convminGetSteps(); + for(Double_t y=convmin;y ¶m){ + +assert(param.size() == 4); + +// check that param are new and hence a calculation is needed (except the phase, which is not needed here) +Bool_t newParams = false; +if (fPrevParam.size() == 0) { + for (UInt_t i=0; i0 +if(width<0) {width=-width;} + +//sort omegas: +if(omega1>=omega2){ + if(omega1>=omega3){ + fOmegaMax=omega1; + if (omega2>=omega3){ + fOmegaCenter=omega2; + fOmegaMin=omega3; + }else{ + fOmegaCenter=omega3; + fOmegaMin=omega2; + } + }else { + fOmegaMin=omega2; + fOmegaCenter=omega1; + fOmegaMax=omega3; + } +}else if (omega2>omega3){ + fOmegaMax=omega2; + if(omega3fOmegaMin){convmax=x-fOmegaMin;} +if(x-convmin0)) fNormalization=1; +return; +} + + +/******************************************************/ +// Implementation PowderLineAxialLor +ClassImp(PowderLineAsymLor) // for the ROOT dictionary + +PowderLineAsymLor::PowderLineAsymLor(){ + fValid = false; + fInvokedGlobal = false; + fGlobalUserFcn = 0; +} + + +PowderLineAsymLor::~PowderLineAsymLor() +{ + if ((fGlobalUserFcn != 0) && fInvokedGlobal) { + delete fGlobalUserFcn; + fGlobalUserFcn = 0; + } +} + + + +void PowderLineAsymLor::SetGlobalPart(vector &globalPart, UInt_t idx) +{ + fIdxGlobal = static_cast(idx); + + if ((Int_t)globalPart.size() <= fIdxGlobal) { // global user function not present, invoke it + fGlobalUserFcn = new PowderLineAsymLorGlobal(); + if (fGlobalUserFcn == 0) { // global user function object couldn't be invoked -> error + fValid = false; + cerr << endl << ">> PowderLineAxialLor::SetGlobalPart(): **ERROR** Couldn't invoke global user function object, sorry ..." << endl; + } else { // global user function object could be invoked -> resize to global user function vector and keep the pointer to the corresponding object + globalPart.resize(fIdxGlobal+1); + globalPart[fIdxGlobal] = dynamic_cast(fGlobalUserFcn); + fValid = true; + fInvokedGlobal = true; + } + } else { // global user function already present hence just retrieve a pointer to it + fValid = true; + fGlobalUserFcn = (PowderLineAsymLorGlobal*)globalPart[fIdxGlobal]; + } +} + + + +Double_t PowderLineAsymLor::operator()(Double_t x, const vector &par) const { + assert(par.size()==4); // make sure the number of parameters handed to the function is correct + Double_t width=par[3];//width of the Lorentzian + + + //ensure width>0 + if(width<0) width=-width; + + //calculate normalization globally + fGlobalUserFcn->UpdateGlobalPart(par); + + // extract Normalization from the global object + Double_t Norm =fGlobalUserFcn->GetNorm(); + Double_t omega_center=fGlobalUserFcn->GetCenter(); + Double_t omega_min=fGlobalUserFcn->GetMin(); + Double_t omega_max=fGlobalUserFcn->GetMax(); + + + //symmetric situation has a delta function: + if(abs(omega_max-omega_min)<0.01*abs(width)){ + cout<omega_min){convmax=x-omega_min;} + if(x-convminGetSteps(); + for(Double_t y=convmin;y ¶m){ + +assert(param.size() == 4); + +// check that param are new and hence a calculation is needed (except the phase, which is not needed here) +Bool_t newParams = false; +if (fPrevParam.size() == 0) { + for (UInt_t i=0; i0 +if(width<0) {width=-width;} + +//sort omegas: +if(omega1>=omega2){ + if(omega1>=omega3){ + fOmegaMax=omega1; + if (omega2>=omega3){ + fOmegaCenter=omega2; + fOmegaMin=omega3; + }else{ + fOmegaCenter=omega3; + fOmegaMin=omega2; + } + }else { + fOmegaMin=omega2; + fOmegaCenter=omega1; + fOmegaMax=omega3; + } +}else if (omega2>omega3){ + fOmegaMax=omega2; + if(omega3fOmegaMin){convmax=x-fOmegaMin;} +if(x-convmin0)) fNormalization=1; +return; +} + + +/******************************************************/ +// Implementation PowderLineAsymGss +ClassImp(PowderLineAsymGss) // for the ROOT dictionary + +PowderLineAsymGss::PowderLineAsymGss(){ + fValid = false; + fInvokedGlobal = false; + fGlobalUserFcn = 0; +} + + +PowderLineAsymGss::~PowderLineAsymGss() +{ + if ((fGlobalUserFcn != 0) && fInvokedGlobal) { + delete fGlobalUserFcn; + fGlobalUserFcn = 0; + } +} + + + +void PowderLineAsymGss::SetGlobalPart(vector &globalPart, UInt_t idx) +{ + fIdxGlobal = static_cast(idx); + + if ((Int_t)globalPart.size() <= fIdxGlobal) { // global user function not present, invoke it + fGlobalUserFcn = new PowderLineAsymGssGlobal(); + if (fGlobalUserFcn == 0) { // global user function object couldn't be invoked -> error + fValid = false; + cerr << endl << ">> PowderLineAxialLor::SetGlobalPart(): **ERROR** Couldn't invoke global user function object, sorry ..." << endl; + } else { // global user function object could be invoked -> resize to global user function vector and keep the pointer to the corresponding object + globalPart.resize(fIdxGlobal+1); + globalPart[fIdxGlobal] = dynamic_cast(fGlobalUserFcn); + fValid = true; + fInvokedGlobal = true; + } + } else { // global user function already present hence just retrieve a pointer to it + fValid = true; + fGlobalUserFcn = (PowderLineAsymGssGlobal*)globalPart[fIdxGlobal]; + } +} + + + +Double_t PowderLineAsymGss::operator()(Double_t x, const vector &par) const { + assert(par.size()==4); // make sure the number of parameters handed to the function is correct + Double_t width=par[3];//width of the Lorentzian + + + //ensure width>0 + if(width<0) width=-width; + + //calculate normalization globally + fGlobalUserFcn->UpdateGlobalPart(par); + + // extract Normalization from the global object + Double_t Norm =fGlobalUserFcn->GetNorm(); + Double_t omega_center=fGlobalUserFcn->GetCenter(); + Double_t omega_min=fGlobalUserFcn->GetMin(); + Double_t omega_max=fGlobalUserFcn->GetMax(); + + + //symmetric situation has a delta function: + if(abs(omega_max-omega_min)<0.01*abs(width)){ + cout<omega_min){convmax=x-omega_min;} + if(x-convminGetSteps(); + for(Double_t y=convmin;y +#include +#include +#include +#include + + +#ifndef PI +#define PI 3.14159265358979323846 /* pi */ +#endif //PI + +#ifndef LINEPROFILE_H +#define LINEPROFILE_H + +using namespace std; + +//helper functions: +Double_t GaussianShape(Double_t, Double_t, Double_t); +Double_t LaplacianShape(Double_t, Double_t, Double_t); +Double_t LorentzianShape(Double_t, Double_t, Double_t); +Double_t IAxial(Double_t, Double_t, Double_t); +Double_t IAsym(Double_t, Double_t, Double_t,Double_t); +Double_t IAsym_low(Double_t, Double_t, Double_t, Double_t); +Double_t IAsym_high(Double_t, Double_t, Double_t, Double_t); + + +// +class LineGauss : public PUserFcnBase { + +public: + + // global user-function-access functions, here without any functionality + Bool_t NeedGlobalPart() const { return false; } + void SetGlobalPart(vector &globalPart, UInt_t idx) { } + Bool_t GlobalPartIsValid() const { return true; } + + // function operator + Double_t operator()(Double_t, const vector&) const; + + // definition of the class for the ROOT dictionary + ClassDef(LineGauss,1) +}; + + +class LineLaplace : public PUserFcnBase { + +public: + + // global user-function-access functions, here without any functionality + Bool_t NeedGlobalPart() const { return false; } + void SetGlobalPart(vector &globalPart, UInt_t idx) { } + Bool_t GlobalPartIsValid() const { return true; } + + // function operator + Double_t operator()(Double_t, const vector&) const; + + // definition of the class for the ROOT dictionary + ClassDef(LineLaplace,1) +}; + + +class LineLorentzian : public PUserFcnBase { + +public: + // global user-function-access functions, here without any functionality + Bool_t NeedGlobalPart() const { return false; } + void SetGlobalPart(vector &globalPart, UInt_t idx) { } + Bool_t GlobalPartIsValid() const { return true; } + + // function operator + Double_t operator()(Double_t, const vector&) const; + + // definition of the class for the ROOT dictionary + ClassDef(LineLorentzian,1) +}; + + +class LineSkewLorentzian : public PUserFcnBase { + +public: + // global user-function-access functions, here without any functionality + Bool_t NeedGlobalPart() const { return false; } + void SetGlobalPart(vector &globalPart, UInt_t idx) { } + Bool_t GlobalPartIsValid() const { return true; } + + // function operator + Double_t operator()(Double_t, const vector&) const; + + // definition of the class for the ROOT dictionary + ClassDef(LineSkewLorentzian,1) +}; + + + +class LineSkewLorentzian2 : public PUserFcnBase { + +public: + // global user-function-access functions, here without any functionality + Bool_t NeedGlobalPart() const { return false; } + void SetGlobalPart(vector &globalPart, UInt_t idx) { } + Bool_t GlobalPartIsValid() const { return true; } + + // function operator + Double_t operator()(Double_t, const vector&) const; + + // definition of the class for the ROOT dictionary + ClassDef(LineSkewLorentzian2,1) +}; + + + +class PowderLineAxialLorGlobal { + +public: + // default constructor and destructor + PowderLineAxialLorGlobal(){} + virtual~PowderLineAxialLorGlobal(){} + + Bool_t IsValid() { return fValid; } + + // the following function will check if something needs to be calculated, which + // is the case if param != fPrevParam + void CalcNormalization(const vector &); + + // this routine returns the globally calculated values + Double_t GetNormalization() const {return fNormalization;} + UInt_t GetSteps() const {return 1e5;} + +private: + Bool_t fValid=true; + vector fPrevParam; + + Double_t fNormalization; + // definition of the class for the ROOT-dictionary + ClassDef(PowderLineAxialLorGlobal,1) +}; + +class PowderLineAxialLor : public PUserFcnBase { + +public: + // default constructor and destructor + PowderLineAxialLor(); + virtual ~PowderLineAxialLor(); + + // global user-function-access functions, here with some functionality + Bool_t NeedGlobalPart() const { return true; } + void SetGlobalPart(vector &globalPart, UInt_t idx); + Bool_t GlobalPartIsValid() const {return fValid;} + + // function operator + Double_t operator()(Double_t, const vector&) const; + +private: + Bool_t fValid; + Bool_t fInvokedGlobal; + Int_t fIdxGlobal; + PowderLineAxialLorGlobal * fGlobalUserFcn; + + // definition of the class for the ROOT dictionary + ClassDef(PowderLineAxialLor,1) +}; + + +class PowderLineAxialGssGlobal { + +public: + // default constructor and destructor + PowderLineAxialGssGlobal(){} + virtual~PowderLineAxialGssGlobal(){} + + Bool_t IsValid() { return fValid; } + + // the following function will check if something needs to be calculated, which + // is the case if param != fPrevParam + void CalcNormalization(const vector &); + + // this routine returns the globally calculated values + Double_t GetNormalization() const {return fNormalization;} + UInt_t GetSteps() const {return 1e5;} + +private: + Bool_t fValid=true; + vector fPrevParam; + + Double_t fNormalization; + // definition of the class for the ROOT-dictionary + ClassDef(PowderLineAxialGssGlobal,1) +}; + +class PowderLineAxialGss : public PUserFcnBase { + +public: + // default constructor and destructor + PowderLineAxialGss(); + virtual ~PowderLineAxialGss(); + + // global user-function-access functions, here with some functionality + Bool_t NeedGlobalPart() const { return true; } + void SetGlobalPart(vector &globalPart, UInt_t idx); + Bool_t GlobalPartIsValid() const {return fValid;} + + // function operator + Double_t operator()(Double_t, const vector&) const; + +private: + Bool_t fValid; + Bool_t fInvokedGlobal; + Int_t fIdxGlobal; + PowderLineAxialGssGlobal * fGlobalUserFcn; + + // definition of the class for the ROOT dictionary + ClassDef(PowderLineAxialGss,1) +}; + + + +class PowderLineAsymLorGlobal { + +public: + // default constructor and destructor + PowderLineAsymLorGlobal(){} + virtual~PowderLineAsymLorGlobal(){} + + Bool_t IsValid() { return fValid; } + + // the following function will check if something needs to be calculated, which + // is the case if param != fPrevParam + void UpdateGlobalPart(const vector &); + + // this routine returns the globally calculated values + Double_t GetNorm() const {return fNormalization;} + Double_t GetCenter() const {return fOmegaCenter;} + Double_t GetMin() const {return fOmegaMin;} + Double_t GetMax() const {return fOmegaMax;} + UInt_t GetSteps() const {return 1e5;} + +private: + Bool_t fValid=true; + vector fPrevParam; + Double_t fNormalization; + Double_t fOmegaCenter; + Double_t fOmegaMin; + Double_t fOmegaMax; + // definition of the class for the ROOT-dictionary + ClassDef(PowderLineAsymLorGlobal,1) +}; + +class PowderLineAsymLor : public PUserFcnBase { + +public: + // default constructor and destructor + PowderLineAsymLor(); + virtual ~PowderLineAsymLor(); + + // global user-function-access functions, here with some functionality + Bool_t NeedGlobalPart() const { return true; } + void SetGlobalPart(vector &globalPart, UInt_t idx); + Bool_t GlobalPartIsValid() const {return fValid;} + + // function operator + Double_t operator()(Double_t, const vector&) const; + +private: + Bool_t fValid; + Bool_t fInvokedGlobal; + Int_t fIdxGlobal; + PowderLineAsymLorGlobal * fGlobalUserFcn; + + // definition of the class for the ROOT dictionary + ClassDef(PowderLineAsymLor,1) +}; + + +class PowderLineAsymGssGlobal { + +public: + // default constructor and destructor + PowderLineAsymGssGlobal(){} + virtual~PowderLineAsymGssGlobal(){} + + Bool_t IsValid() { return fValid; } + + // the following function will check if something needs to be calculated, which + // is the case if param != fPrevParam + void UpdateGlobalPart(const vector &); + + // this routine returns the globally calculated values + Double_t GetNorm() const {return fNormalization;} + Double_t GetCenter() const {return fOmegaCenter;} + Double_t GetMin() const {return fOmegaMin;} + Double_t GetMax() const {return fOmegaMax;} + UInt_t GetSteps() const {return 1e5;} + +private: + Bool_t fValid=true; + vector fPrevParam; + Double_t fNormalization; + Double_t fOmegaCenter; + Double_t fOmegaMin; + Double_t fOmegaMax; + // definition of the class for the ROOT-dictionary + ClassDef(PowderLineAsymGssGlobal,1) +}; + +class PowderLineAsymGss : public PUserFcnBase { + +public: + // default constructor and destructor + PowderLineAsymGss(); + virtual ~PowderLineAsymGss(); + + // global user-function-access functions, here with some functionality + Bool_t NeedGlobalPart() const { return true; } + void SetGlobalPart(vector &globalPart, UInt_t idx); + Bool_t GlobalPartIsValid() const {return fValid;} + + // function operator + Double_t operator()(Double_t, const vector&) const; + +private: + Bool_t fValid; + Bool_t fInvokedGlobal; + Int_t fIdxGlobal; + PowderLineAsymGssGlobal * fGlobalUserFcn; + + // definition of the class for the ROOT dictionary + ClassDef(PowderLineAsymGss,1) +}; + + + + +#endif //LINEPROFILE_H diff --git a/src/external/libBNMR/libLineProfile/libLineProfileLinkDef.h b/src/external/libBNMR/libLineProfile/libLineProfileLinkDef.h new file mode 100644 index 00000000..50d55a53 --- /dev/null +++ b/src/external/libBNMR/libLineProfile/libLineProfileLinkDef.h @@ -0,0 +1,50 @@ +/**************************************************************************/ + +/*************************************************************************** + * Copyright (C) 2017 by Jonas A. Krieger * + * jonas.krieger@psi.ch * + * * + * This program is free software; you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation; either version 2 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program; if not, write to the * + * Free Software Foundation, Inc., * + * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. * +***************************************************************************/ + +//root dictionary +#ifdef __CINT__ + +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; + +#pragma link C++ class LineGauss+; +#pragma link C++ class LineLaplace+; +#pragma link C++ class LineLorentzian+; +#pragma link C++ class LineSkewLorentzian+; +#pragma link C++ class LineSkewLorentzian2+; +#pragma link C++ class PowderLineAxialLor+; +#pragma link C++ class PowderLineAxialLorGlobal+; +#pragma link C++ class PowderLineAxialGss+; +#pragma link C++ class PowderLineAxialGssGlobal+; +#pragma link C++ class PowderLineAsymLor+; +#pragma link C++ class PowderLineAsymLorGlobal+; +#pragma link C++ class PowderLineAsymGss+; +#pragma link C++ class PowderLineAsymGssGlobal+; + + + + + + + +#endif //__CINT__ diff --git a/src/external/libPhotoMeissner/classes/CMakeLists.txt b/src/external/libPhotoMeissner/classes/CMakeLists.txt index e52f90b0..dcf94dec 100644 --- a/src/external/libPhotoMeissner/classes/CMakeLists.txt +++ b/src/external/libPhotoMeissner/classes/CMakeLists.txt @@ -4,8 +4,10 @@ set(MUSRFIT_INC ${CMAKE_SOURCE_DIR}/src/include) root_generate_dictionary( - PPhotoMeissnerDict + PPhotoMeissnerDict -I${FFTW3_INCLUDE_DIR} + -I${GSL_INCLUDE_DIRS} + -I${ROOT_INCLUDE_DIRS} -I${MUSRFIT_INC} -I${CMAKE_CURRENT_SOURCE_DIR}/../include PPhotoMeissner.h @@ -46,8 +48,10 @@ set_target_properties(PPhotoMeissner #--- make sure that the include directory is found ---------------------------- target_include_directories( - PPhotoMeissner BEFORE PRIVATE + PPhotoMeissner BEFORE PRIVATE $ + $ + $ $ $ )