Added functions for beta-NMR relaxation curves (beam on/off).

This commit is contained in:
salman 2010-09-02 18:50:11 +00:00
parent 6774573f7a
commit cd9a625958
5 changed files with 286 additions and 0 deletions

60
src/external/libBNMR/Makefile.libBNMR vendored Normal file
View File

@ -0,0 +1,60 @@
#---------------------------------------------------
# get compilation flags from root-config
ROOTCFLAGS = $(shell $(ROOTSYS)/bin/root-config --cflags)
#---------------------------------------------------
OS = LINUX
CXX = g++
CXXFLAGS = -O3 -Wall -Wno-trigraphs -fPIC
LOCALINCLUDE = .
ROOTINCLUDE = $(ROOTSYS)/include
INCLUDES = -I$(LOCALINCLUDE) -I$(ROOTINCLUDE)
LD = g++
LDFLAGS =
SOFLAGS = -O -shared
# the output from the root-config script:
CXXFLAGS += $(ROOTCFLAGS)
LDFLAGS +=
# some definitions: headers (used to generate *Dict* stuff), sources, objects,...
OBJS =
OBJS += TBNMR.o TlibBNMRDict.o
SHLIB = libBNMR.so
# make the shared lib:
#
all: $(SHLIB)
$(SHLIB): $(OBJS)
@echo "---> 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
TlibBNMRDict.cpp: TBNMR.h TBNMRLinkDef.h
@echo "Generating dictionary $@..."
rootcint -f $@ -c -p -I$(ROOTINCLUDE) $^
install: all
@echo "Installing shared lib: libTApproximation.so"
ifeq ($(OS),LINUX)
cp -pv $(SHLIB) $(ROOTSYS)/lib
cp -pv $(LOCALINCLUDE)/*.h $(ROOTSYS)/include
endif

104
src/external/libBNMR/TBNMR.cpp vendored Normal file
View File

@ -0,0 +1,104 @@
/***************************************************************************
TBNMR.cpp
Author: Zaher Salman
e-mail: zaher.salman@psi.ch
2010/09/02
***************************************************************************/
/***************************************************************************
* Copyright (C) 2010 by Zaher Salman *
* zaher.salman@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 "TBNMR.h"
#include "TF1.h"
#include "Math/WrappedTF1.h"
#include "Math/GaussIntegrator.h"
#define tau_Li 1210
#define gamma_Li 6.3018 // In units kHz/mT
#define PI 3.14159265358979323846
#define TWOPI 6.28318530717958647692
ClassImp(TBNMR) // for the ROOT-dictionary
ClassImp(ExpRlx)
ClassImp(SExpRlx)
double TBNMR::operator()(double x, const vector<double> &par) const {
assert(par.size()==1); // make sure the number of parameters handed to the function is correct
double arg(par[0]*x);
if(!arg)
return 1.0;
return sin(arg)/arg;
}
double ExpRlx::operator()(double x, const vector<double> &par) const {
assert(par.size()==2); // make sure the number of parameters handed to the function is correct
// par[0] time of beam off
// par[1] is the relaxation rate
double tau_p;
double y;
tau_p = (tau_Li/(1.+par[1]*tau_Li));
if ( x <= par[0] && x >= 0) {
y=(tau_p/tau_Li)*(1-exp(-x/tau_p))/(1-exp(-x/tau_Li));
} else if ( x > par[0] ){
y=(tau_p/tau_Li)*(1-exp(-par[0]/tau_p))/(1-exp(-par[0]/tau_Li))*exp(-par[1]*(x-par[0]));
} else {
y = 0;
}
return y;
}
double SExpRlx::operator()(double x, const vector<double> &par) const {
assert(par.size()==3); // make sure the number of parameters handed to the function is correct
// par[0] time of beam off
// par[1] is the relaxation rate
// par[2] is the exponent
double tau_p;
double y;
tau_p = (tau_Li/(1.+par[1]*tau_Li));
if ( x >= 0 && x <= par[0] ) {
TF1 sexp("sexp", "exp(-([0]-x)/[3])*exp(-pow(([1]*([0]-x)),[2]))", 0.0, 10000.0);
sexp.SetParameters(x, par[1], par[2],tau_Li);
sexp.SetNpx(1000);
y=sexp.Integral(0.0,x)/(1-exp(-x/tau_Li))/tau_Li;
} else if ( x > par[0] ) {
TF1 sexp("sexp", "exp(-([3]-x)/[4])*exp(-pow(([1]*([0]-x)),[2]))", 0.0, 10000.0);
sexp.SetParameters(x, par[1], par[2], par[0],tau_Li);
sexp.SetNpx(1000);
y=sexp.Integral(0.0,par[0])/(1-exp(-x/tau_Li))/tau_Li;
} else {
y = 0;
}
return y;
}

80
src/external/libBNMR/TBNMR.h vendored Normal file
View File

@ -0,0 +1,80 @@
/***************************************************************************
TBNMR.h
Author: Zaher Salman
e-mail: zaher.salman@psi.ch
2010/09/02
***************************************************************************/
/***************************************************************************
* Copyright (C) 2010 by Zaher Salman *
* zaher.salman@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 "PUserFcnBase.h"
#include <cassert>
#include <cmath>
#include <vector>
using namespace std;
class TBNMR : public PUserFcnBase {
public:
// default constructor and destructor
TBNMR(){}
~TBNMR(){}
// function operator
double operator()(double, const vector<double>&) const;
// definition of the class for the ROOT-dictionary
ClassDef(TBNMR,1)
};
class ExpRlx : public PUserFcnBase {
public:
// default constructor and destructor
ExpRlx(){}
~ExpRlx(){}
// function operator
double operator()(double, const vector<double>&) const;
// definition of the class for the ROOT-dictionary
ClassDef(ExpRlx,1)
};
class SExpRlx : public PUserFcnBase {
public:
// default constructor and destructor
SExpRlx(){}
~SExpRlx(){}
// function operator
double operator()(double, const vector<double>&) const;
// definition of the class for the ROOT-dictionary
ClassDef(SExpRlx,1)
};

BIN
src/external/libBNMR/TBNMR.o vendored Normal file

Binary file not shown.

42
src/external/libBNMR/TBNMRLinkDef.h vendored Normal file
View File

@ -0,0 +1,42 @@
/***************************************************************************
TBNMRLinkDef.h
Author: Zaher Salman
e-mail: zaher.salman@psi.ch
2010/09/02
***************************************************************************/
/***************************************************************************
* Copyright (C) 2010 by Zaher Salman *
* zaher.salman@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. *
***************************************************************************/
#ifdef __CINT__
#pragma link off all globals;
#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+;
#endif //__CINT__