LMU/musrfit
LMU/musrfit
5
0
Fork 0
Code Issues 6 Pull Requests Actions Packages Projects Releases 3 Activity

changed the user function for nonlocal fitting to the new rge-handler. This breaks the backwards compatibility!

Browse Source
This commit is contained in:
suter_a 2021-03-30 15:48:51 +02:00
parent df53338d6b
commit a821890ef5
11 changed files with 90 additions and 474 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
doc/examples/ASlibs/nonlocal_startup.xml
View File

@ -1,13 +1,14 @@
<?xml version="1.0" encoding="UTF-8"?> <?xml version="1.0" encoding="UTF-8"?>
<nonlocal xmlns="http://nemu.web.psi.ch/musrfit/nonlocal"> <nonlocal xmlns="http://nemu.web.psi.ch/musrfit/nonlocal">
<comment> <comment>
$Id: nonlocal_startup.xml 4047 2009-07-02 13:36:18Z nemu $ nonlocal_startup.xml
</comment> </comment>
<nonlocal_par> <nonlocal_par>
<fourier_points>262144</fourier_points> <fourier_points>262144</fourier_points>
</nonlocal_par> </nonlocal_par>
<trim_sp_part> <trim_sp>
<data_path>profiles/Sn_E</data_path> <data_path>./profiles/</data_path>
<rge_fln_pre>Sn_E</rge_fln_pre>
<energy_list> <energy_list>
<energy>1000</energy> <energy>1000</energy>
<energy>2000</energy> <energy>2000</energy>
@ -22,5 +23,5 @@
<energy>25000</energy> <energy>25000</energy>
<energy>27300</energy> <energy>27300</energy>
</energy_list> </energy_list>
</trim_sp_part> </trim_sp>
</nonlocal> </nonlocal>

32
doc/examples/ASlibs/test-nonlocal.msr
View File

@ -4,21 +4,21 @@ FITPARAMETER
# Nr. Name Value Step Pos_Error Boundaries # Nr. Name Value Step Pos_Error Boundaries
1 one 1 0 none 1 one 1 0 none
2 zero 0 0 none 2 zero 0 0 none
3 Asy 0.2235 -0.0066 0.0072 0 0.33 3 Asy 0.2235 -0.0065 0.0072 0 0.33
4 energy 22 0 none 4 energy 22000 0 none
5 redTemp 0.8683 0 none 5 redTemp 0.8683 0 none
6 thickness 5000 0 none 6 thickness 5000 0 none
7 ell 12000 0 none 7 ell 12000 0 none
8 xi 94 0 none 8 xi 94 0 none
9 lambdaL 52.5 -1.4 1.2 9 lambdaL 52.5 -1.3 1.2
10 Bext 47.11 0 none 10 Bext 47.11 0 none
11 deadLayer 0.00052 -0.00052 6.89563 0 none 11 deadLayer 0.0004 0.0037 none 0 none
12 RateSmear 0.418 -0.015 0.015 12 RateSmear 0.418 -0.015 0.015
13 N0_L 301.61 0.62 none 13 N0_L 301.6 -1.1 1.2
14 Bkg_L 24.442 0.076 none 14 Bkg_L 24.441 -0.092 0.092
15 Phase_L 42.0 -2.2 2.2 15 Phase_L 42.0 -2.2 2.2
16 Alpha_LR 1.0614 -0.0066 0.0064 16 Alpha_LR 1.0613 -0.0066 0.0064
17 Bkg_R 27.646 -0.090 0.090 17 Bkg_R 27.646 -0.089 0.090
18 RelPhase_R 154.4 -1.1 1.1 18 RelPhase_R 154.4 -1.1 1.1
############################################################### ###############################################################
@ -37,7 +37,6 @@ RUN data/lem10_his_0825 MUE4 PSI ROOT-NPP (name beamline institute data-file-f
fittype 0 (single histogram fit) fittype 0 (single histogram fit)
norm 13 norm 13
backgr.fit 14 backgr.fit 14
lifetimecorrection
map 1 2 0 0 0 0 0 0 0 0 map 1 2 0 0 0 0 0 0 0 0
forward 1 forward 1
data 3289 65000 data 3289 65000
@ -49,7 +48,6 @@ RUN data/lem10_his_0825 MUE4 PSI ROOT-NPP (name beamline institute data-file-f
fittype 0 (single histogram fit) fittype 0 (single histogram fit)
norm fun1 norm fun1
backgr.fit 17 backgr.fit 17
lifetimecorrection
map 16 18 0 0 0 0 0 0 0 0 map 16 18 0 0 0 0 0 0 0 0
forward 3 forward 3
data 3289 65000 data 3289 65000
@ -69,23 +67,27 @@ SAVE
############################################################### ###############################################################
PLOT 0 (single histo plot) PLOT 0 (single histo plot)
lifetimecorrection
runs 1 2 runs 1 2
range 0 9 -0.3 0.3 range 0 9 -0.3 0.3
view_packing 500 view_packing 500
############################################################### ###############################################################
FOURIER FOURIER
units Gauss # units either 'Gauss', 'MHz', or 'Mc/s' units Gauss # units either 'Gauss', 'Tesla', 'MHz', or 'Mc/s'
fourier_power 10 fourier_power 10
apodization STRONG # NONE, WEAK, MEDIUM, STRONG apodization STRONG # NONE, WEAK, MEDIUM, STRONG
plot POWER # REAL, IMAG, REAL_AND_IMAG, POWER, PHASE plot POWER # REAL, IMAG, REAL_AND_IMAG, POWER, PHASE, PHASE_OPT_REAL
phase 8.5 phase 8.5
#range_for_phase_correction 50.0 70.0 #range_for_phase_correction 50.0 70.0
range 0.0 200.0 range 0 200
############################################################### ###############################################################
STATISTIC --- 2014-11-05 12:34:54 STATISTIC --- 2021-03-30 15:04:30
maxLH = 415.6, NDF = 358, maxLH/NDF = 1.160933 maxLH = 415.6, NDF = 358, maxLH/NDF = 1.160942
expected maxLH = 418.1, NDF = 358, expected maxLH/NDF = 1.168011
run block 1: (NDF/red.maxLH/red.maxLH_e) = (177/1.072511/1.087460)
run block 2: (NDF/red.maxLH/red.maxLH_e) = (175/1.290186/1.289528)

23
src/external/Nonlocal/CMakeLists.txt vendored
View File

@ -2,6 +2,7 @@
#--- generate necessary dictionaries ------------------------------------------ #--- generate necessary dictionaries ------------------------------------------
set(MUSRFIT_INC ${CMAKE_SOURCE_DIR}/src/include) set(MUSRFIT_INC ${CMAKE_SOURCE_DIR}/src/include)
set(MUSRFIT_CLASSES ${CMAKE_SOURCE_DIR}/src/classes)
set(NONLOCAL_INC ${CMAKE_SOURCE_DIR}/src/external/Nonlocal) set(NONLOCAL_INC ${CMAKE_SOURCE_DIR}/src/external/Nonlocal)
# ROOT requires that the dictonary header files are found at configuration time. # ROOT requires that the dictonary header files are found at configuration time.
# Hence, target_include_directories cannot be used here because, targets are # Hence, target_include_directories cannot be used here because, targets are
@ -30,13 +31,24 @@ root_generate_dictionary(
LINKDEF PNL_StartupHandlerLinkDef.h LINKDEF PNL_StartupHandlerLinkDef.h
MODULE PNL_StartupHandler MODULE PNL_StartupHandler
) )
root_generate_dictionary(
PRgeHandlerDict
PRgeHandler.h
OPTIONS
-I${NONLOCAL_INC}
-I${MUSRFIT_INC}
-I${CMAKE_CURRENT_SOURCE_DIR}
-inlineInputHeader
LINKDEF ${MUSRFIT_INC}/PRgeHandlerLinkDef.h
MODULE PRgeHandler
)
#--- create pkg-config info --------------------------------------------------- #--- create pkg-config info ---------------------------------------------------
set(prefix "${CMAKE_INSTALL_PREFIX}") set(prefix "${CMAKE_INSTALL_PREFIX}")
set(exec_prefix "\$\{prefix\}") set(exec_prefix "\$\{prefix\}")
set(libdir "\$\{exec_prefix\}/lib") set(libdir "\$\{exec_prefix\}/lib")
set(includedir "\$\{prefix\}/include") set(includedir "\$\{prefix\}/include")
set(PNL_PIPPARDFITTER_VERSION "1.0.0") set(PNL_PIPPARDFITTER_VERSION "1.1.0")
set(PNL_PIPPARDFITTER_LIBRARY_NAME "PNL_PippardFitter") set(PNL_PIPPARDFITTER_LIBRARY_NAME "PNL_PippardFitter")
configure_file("PNL_PippardFitter.pc.in" "PNL_PippardFitter.pc" @ONLY) configure_file("PNL_PippardFitter.pc.in" "PNL_PippardFitter.pc" @ONLY)
@ -44,7 +56,8 @@ configure_file("PNL_PippardFitter.pc.in" "PNL_PippardFitter.pc" @ONLY)
add_library(PNL_PippardFitter SHARED add_library(PNL_PippardFitter SHARED
PNL_PippardFitter.cpp PNL_PippardFitter.cpp
PNL_PippardFitterDict.cxx PNL_PippardFitterDict.cxx
PNL_RgeHandler.cpp ${MUSRFIT_CLASSES}/PRgeHandler.cpp
PRgeHandlerDict.cxx
PNL_StartupHandler.cpp PNL_StartupHandler.cpp
PNL_StartupHandlerDict.cxx PNL_StartupHandlerDict.cxx
) )
@ -64,7 +77,7 @@ set_target_properties(PNL_PippardFitter
) )
#--- add library dependencies ------------------------------------------------- #--- add library dependencies -------------------------------------------------
target_link_libraries(PNL_PippardFitter ${FFTW3_LIBRARY} ${ROOT_LIBRARIES} PUserFcnBase) target_link_libraries(PNL_PippardFitter ${Boost_LIBRARIES} ${FFTW3_LIBRARY} ${ROOT_LIBRARIES} PUserFcnBase)
#--- install PNL_PippardFitter solib ------------------------------------------ #--- install PNL_PippardFitter solib ------------------------------------------
install(TARGETS PNL_PippardFitter DESTINATION lib) install(TARGETS PNL_PippardFitter DESTINATION lib)
@ -74,6 +87,8 @@ install(
FILES ${CMAKE_CURRENT_BINARY_DIR}/libPNL_PippardFitter_rdict.pcm FILES ${CMAKE_CURRENT_BINARY_DIR}/libPNL_PippardFitter_rdict.pcm
${CMAKE_CURRENT_BINARY_DIR}/libPNL_PippardFitter.rootmap ${CMAKE_CURRENT_BINARY_DIR}/libPNL_PippardFitter.rootmap
${CMAKE_CURRENT_BINARY_DIR}/libPNL_StartupHandler_rdict.pcm ${CMAKE_CURRENT_BINARY_DIR}/libPNL_StartupHandler_rdict.pcm
${CMAKE_CURRENT_BINARY_DIR}/libPRgeHandler.rootmap
${CMAKE_CURRENT_BINARY_DIR}/libPRgeHandler_rdict.pcm
${CMAKE_CURRENT_BINARY_DIR}/libPNL_StartupHandler.rootmap ${CMAKE_CURRENT_BINARY_DIR}/libPNL_StartupHandler.rootmap
DESTINATION lib DESTINATION lib
) )
@ -81,9 +96,7 @@ install(
#--- install PNL_PippardFitter header ----------------------------------------- #--- install PNL_PippardFitter header -----------------------------------------
install( install(
FILES FILES
PNonlocal.h
PNL_PippardFitter.h PNL_PippardFitter.h
PNL_RgeHandler.h
PNL_StartupHandler.h PNL_StartupHandler.h
DESTINATION DESTINATION
include include

35
src/external/Nonlocal/PNL_PippardFitter.cpp vendored
View File

@ -31,7 +31,6 @@
#include <cmath> #include <cmath>
#include <iostream> #include <iostream>
using namespace std;
#include <TSAXParser.h> #include <TSAXParser.h>
#include <TMath.h> #include <TMath.h>
@ -67,8 +66,8 @@ PNL_PippardFitterGlobal::PNL_PippardFitterGlobal()
Int_t status = parseXmlFile(saxParser, startup_path_name); Int_t status = parseXmlFile(saxParser, startup_path_name);
// check for parse errors // check for parse errors
if (status) { // error if (status) { // error
cout << endl << ">> PNL_PippardFitterGlobal::PNL_PippardFitterGlobal: **WARNING** Reading/parsing nonlocal_startup.xml failed."; std::cout << std::endl << ">> PNL_PippardFitterGlobal::PNL_PippardFitterGlobal: **WARNING** Reading/parsing nonlocal_startup.xml failed.";
cout << endl; std::cout << std::endl;
fValid = false; fValid = false;
} }
@ -80,20 +79,20 @@ PNL_PippardFitterGlobal::PNL_PippardFitterGlobal()
// check if everything went fine with the startup handler // check if everything went fine with the startup handler
if (!fStartupHandler->IsValid()) { if (!fStartupHandler->IsValid()) {
cout << endl << ">> PNL_PippardFitterGlobal::PNL_PippardFitterGlobal **PANIC ERROR**"; std::cout << std::endl << ">> PNL_PippardFitterGlobal::PNL_PippardFitterGlobal **PANIC ERROR**";
cout << endl << ">> startup handler too unhappy. Will terminate unfriendly, sorry."; std::cout << std::endl << ">> startup handler too unhappy. Will terminate unfriendly, sorry.";
cout << endl; std::cout << std::endl;
fValid = false; fValid = false;
} }
fFourierPoints = fStartupHandler->GetFourierPoints(); fFourierPoints = fStartupHandler->GetFourierPoints();
// load all the TRIM.SP rge-files // load all the TRIM.SP rge-files
fRgeHandler = new PNL_RgeHandler(fStartupHandler->GetTrimSpDataPathList(), fStartupHandler->GetTrimSpDataVectorList()); fRgeHandler = new PRgeHandler("./nonlocal_startup.xml");
if (!fRgeHandler->IsValid()) { if (!fRgeHandler->IsValid()) {
cout << endl << ">> PNL_PippardFitterGlobal::PNL_PippardFitterGlobal **PANIC ERROR**"; std::cout << std::endl << ">> PNL_PippardFitterGlobal::PNL_PippardFitterGlobal **PANIC ERROR**";
cout << endl << ">> rge data handler too unhappy. Will terminate unfriendly, sorry."; std::cout << std::endl << ">> rge data handler too unhappy. Will terminate unfriendly, sorry.";
cout << endl; std::cout << std::endl;
fValid = false; fValid = false;
} }
@ -182,16 +181,16 @@ void PNL_PippardFitterGlobal::CalculateField(const std::vector<Double_t> &param)
if (fFieldq == 0) { if (fFieldq == 0) {
fFieldq = (fftw_complex *) fftw_malloc(sizeof(fftw_complex) * fFourierPoints); fFieldq = (fftw_complex *) fftw_malloc(sizeof(fftw_complex) * fFourierPoints);
if (fFieldq == 0) { if (fFieldq == 0) {
cout << endl << "PPippard::CalculateField(): **ERROR** couldn't allocate memory for fFieldq"; std::cout << std::endl << "PPippard::CalculateField(): **ERROR** couldn't allocate memory for fFieldq";
cout << endl; std::cout << std::endl;
return; return;
} }
} }
if (fFieldB == 0) { if (fFieldB == 0) {
fFieldB = (fftw_complex *) fftw_malloc(sizeof(fftw_complex) * fFourierPoints); fFieldB = (fftw_complex *) fftw_malloc(sizeof(fftw_complex) * fFourierPoints);
if (fFieldB == 0) { if (fFieldB == 0) {
cout << endl << "PPippard::CalculateField(): **ERROR** couldn't allocate memory for fFieldB"; std::cout << std::endl << "PPippard::CalculateField(): **ERROR** couldn't allocate memory for fFieldB";
cout << endl; std::cout << std::endl;
return; return;
} }
} }
@ -387,7 +386,7 @@ PNL_PippardFitter::~PNL_PippardFitter()
* \param globalPart * \param globalPart
* \param idx * \param idx
*/ */
void PNL_PippardFitter::SetGlobalPart(vector<void*> &globalPart, UInt_t idx) void PNL_PippardFitter::SetGlobalPart(std::vector<void*> &globalPart, UInt_t idx)
{ {
fIdxGlobal = static_cast<Int_t>(idx); fIdxGlobal = static_cast<Int_t>(idx);
@ -395,10 +394,10 @@ void PNL_PippardFitter::SetGlobalPart(vector<void*> &globalPart, UInt_t idx)
fPippardFitterGlobal = new PNL_PippardFitterGlobal(); fPippardFitterGlobal = new PNL_PippardFitterGlobal();
if (fPippardFitterGlobal == 0) { if (fPippardFitterGlobal == 0) {
fValid = false; fValid = false;
cerr << endl << ">> PNL_PippardFitter::SetGlobalPart(): **ERROR** Couldn't invoke global user function object, sorry ..." << endl; std::cerr << std::endl << ">> PNL_PippardFitter::SetGlobalPart(): **ERROR** Couldn't invoke global user function object, sorry ..." << std::endl;
} else if (!fPippardFitterGlobal->IsValid()) { } else if (!fPippardFitterGlobal->IsValid()) {
fValid = false; fValid = false;
cerr << endl << ">> PNL_PippardFitter::SetGlobalPart(): **ERROR** initialization of global user function object failed, sorry ..." << endl; std::cerr << std::endl << ">> PNL_PippardFitter::SetGlobalPart(): **ERROR** initialization of global user function object failed, sorry ..." << std::endl;
} else { } else {
fValid = true; fValid = true;
fInvokedGlobal = true; fInvokedGlobal = true;
@ -447,7 +446,7 @@ Double_t PNL_PippardFitter::operator()(Double_t t, const std::vector<Double_t> &
fPippardFitterGlobal->CalculateField(param); fPippardFitterGlobal->CalculateField(param);
Int_t energyIndex = fPippardFitterGlobal->GetEnergyIndex(param[0]); Int_t energyIndex = fPippardFitterGlobal->GetEnergyIndex(param[0]);
if (energyIndex == -1) { // energy not found if (energyIndex == -1) { // energy not found
cerr << endl << ">> PNL_PippardFitter::operator() energy " << param[0] << " not found. Will terminate." << endl; std::cerr << std::endl << ">> PNL_PippardFitter::operator() energy " << param[0] << " not found. Will terminate." << std::endl;
assert(0); assert(0);
} }

11
src/external/Nonlocal/PNL_PippardFitter.h vendored
View File

@ -36,9 +36,10 @@
//#endif //#endif
#include <fftw3.h> #include <fftw3.h>
#include "PMusr.h"
#include "PUserFcnBase.h" #include "PUserFcnBase.h"
#include "PNL_StartupHandler.h" #include "PNL_StartupHandler.h"
#include "PNL_RgeHandler.h" #include "PRgeHandler.h"
class PNL_PippardFitterGlobal class PNL_PippardFitterGlobal
{ {
@ -49,15 +50,15 @@ class PNL_PippardFitterGlobal
Bool_t IsValid() { return fValid; } Bool_t IsValid() { return fValid; }
virtual void SetTempExponent(const Double_t nn) { f_nn = nn; } virtual void SetTempExponent(const Double_t nn) { f_nn = nn; }
virtual void CalculateField(const std::vector<Double_t> &param) const; virtual void CalculateField(const std::vector<Double_t> &param) const;
virtual Int_t GetEnergyIndex(const Double_t energy) { return fRgeHandler->GetRgeEnergyIndex(energy); } virtual Int_t GetEnergyIndex(const Double_t energy) { return fRgeHandler->GetEnergyIndex(energy); }
virtual Double_t GetMuonStoppingDensity(const Int_t energyIndex, const Double_t z) const { return fRgeHandler->GetRgeValue(energyIndex, z); } virtual Double_t GetMuonStoppingDensity(const Int_t energyIndex, const Double_t z) const { return fRgeHandler->Get_n(energyIndex, z); }
virtual Double_t GetMagneticField(const Double_t z) const; virtual Double_t GetMagneticField(const Double_t z) const;
private: private:
Bool_t fValid; Bool_t fValid;
PNL_StartupHandler *fStartupHandler; PNL_StartupHandler *fStartupHandler;
PNL_RgeHandler *fRgeHandler; PRgeHandler *fRgeHandler;
mutable std::vector<Double_t> fPreviousParam; mutable std::vector<Double_t> fPreviousParam;
@ -89,7 +90,7 @@ class PNL_PippardFitter : public PUserFcnBase
virtual ~PNL_PippardFitter(); virtual ~PNL_PippardFitter();
virtual Bool_t NeedGlobalPart() const { return true; } virtual Bool_t NeedGlobalPart() const { return true; }
virtual void SetGlobalPart(vector<void*> &globalPart, UInt_t idx); virtual void SetGlobalPart(std::vector<void*> &globalPart, UInt_t idx);
virtual Bool_t GlobalPartIsValid() const; virtual Bool_t GlobalPartIsValid() const;
virtual Double_t operator()(Double_t t, const std::vector<Double_t> &param) const; virtual Double_t operator()(Double_t t, const std::vector<Double_t> &param) const;

208
src/external/Nonlocal/PNL_RgeHandler.cpp vendored
View File

@ -1,208 +0,0 @@
/***************************************************************************
PNL_RgeHandler.cpp
Author: Andreas Suter
e-mail: andreas.suter@psi.ch
$Id$
***************************************************************************/
/***************************************************************************
* Copyright (C) 2009 by Andreas Suter *
* andreas.suter@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 <cassert>
#include <iostream>
#include <fstream>
#include "PNL_RgeHandler.h"
//--------------------------------------------------------------------------
// Constructor
//--------------------------------------------------------------------------
/**
*
*/
PNL_RgeHandler::PNL_RgeHandler(const PStringVector &rgeDataPathList, const PDoubleVector &rgeDataEnergyList)
{
fIsValid = false;
fIsValid = LoadRgeData(rgeDataPathList, rgeDataEnergyList);
}
//--------------------------------------------------------------------------
// Destructor
//--------------------------------------------------------------------------
/**
*
*/
PNL_RgeHandler::~PNL_RgeHandler()
{
fRgeDataList.clear();
}
//--------------------------------------------------------------------------
// GetRgeEnergyIndex
//--------------------------------------------------------------------------
/**
*
*/
Int_t PNL_RgeHandler::GetRgeEnergyIndex(const Double_t energy)
{
Int_t idx = -1;
for (UInt_t i=0; i<fRgeDataList.size(); i++) {
if (energy == fRgeDataList[i].energy) {
idx = i;
break;
}
}
return idx;
}
//--------------------------------------------------------------------------
// GetRgeValue
//--------------------------------------------------------------------------
/**
*
*/
Double_t PNL_RgeHandler::GetRgeValue(const Int_t index, const Double_t dist)
{
Double_t rgeVal = 0.0;
// find the proper position index
Int_t distIdx = -2;
for (UInt_t i=0; i<fRgeDataList[index].stoppingDistance.size(); i++) {
if (fRgeDataList[index].stoppingDistance[i] > dist) {
distIdx = i-1;
break;
}
}
// get the proper n(z) value
if (distIdx < 0) {
rgeVal = 0.0;
} else {
rgeVal = fRgeDataList[index].stoppingAmplitude[distIdx] +
(fRgeDataList[index].stoppingAmplitude[distIdx+1] - fRgeDataList[index].stoppingAmplitude[distIdx]) *
(dist-fRgeDataList[index].stoppingDistance[distIdx])/(fRgeDataList[index].stoppingDistance[distIdx+1]-fRgeDataList[index].stoppingDistance[distIdx]);
}
return rgeVal;
}
//--------------------------------------------------------------------------
// GetRgeValue
//--------------------------------------------------------------------------
/**
*
*/
Double_t PNL_RgeHandler::GetRgeValue(const Double_t energy, const Double_t dist)
{
// check if energy is present in rge data list
Int_t idx = -1;
for (UInt_t i=0; i<fRgeDataList.size(); i++) {
if (energy == fRgeDataList[i].energy) {
idx = i;
break;
}
}
// energy not found
if (idx == -1)
return -1.0;
return GetRgeValue(idx, dist);
}
//--------------------------------------------------------------------------
// LoadRgeData
//--------------------------------------------------------------------------
/**
*
*/
Bool_t PNL_RgeHandler::LoadRgeData(const PStringVector &rgeDataPathList, const PDoubleVector &rgeDataEnergyList)
{
std::ifstream fin;
PNL_RgeData data;
TString tstr;
char line[512];
int result;
double dist, val;
for (UInt_t i=0; i<rgeDataPathList.size(); i++) {
// open rge-file for reading
fin.open(rgeDataPathList[i].Data(), std::iostream::in);
if (!fin.is_open()) {
std::cout << std::endl << "PNL_RgeHandler::LoadRgeData **ERROR**";
std::cout << std::endl << " Could not open file " << rgeDataPathList[i].Data();
std::cout << std::endl;
return false;
}
// keep energy (in keV)
data.energy = rgeDataEnergyList[i]/1000.0;
// read msr-file
while (!fin.eof()) {
// read a line
fin.getline(line, sizeof(line));
// ignore empty lines
tstr = line;
if (tstr.IsWhitespace())
continue;
// get values
result = sscanf(line, "%lf %lf", &dist, &val);
// check if data are valid, otherwise ignore it
if (result != 2) {
continue;
}
// feed fRgeDataList
data.stoppingDistance.push_back(dist/10.0);
data.stoppingAmplitude.push_back(val);
}
// normalize stopping distribution
Double_t norm = 0.0;
for (UInt_t j=0; j<data.stoppingAmplitude.size(); j++)
norm += data.stoppingAmplitude[j];
norm *= (data.stoppingDistance[1] - data.stoppingDistance[0]);
for (UInt_t j=0; j<data.stoppingAmplitude.size(); j++)
data.stoppingAmplitude[j] /= norm;
// keep data
fRgeDataList.push_back(data);
// clean up
data.stoppingAmplitude.clear();
data.stoppingDistance.clear();
fin.close();
}
return true;
}

55
src/external/Nonlocal/PNL_RgeHandler.h vendored
View File

@ -1,55 +0,0 @@
/***************************************************************************
PNL_RgeHandler.h
Author: Andreas Suter
e-mail: andreas.suter@psi.ch
$Id$
***************************************************************************/
/***************************************************************************
* Copyright (C) 2009-2021 by Andreas Suter *
* andreas.suter@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. *
***************************************************************************/
#ifndef _PNL_RGEHANDLER_H_
#define _PNL_RGEHANDLER_H_
#include "PNonlocal.h"
class PNL_RgeHandler
{
public:
PNL_RgeHandler(const PStringVector &rgeDataPathList, const PDoubleVector &rgeDataEnergyList);
virtual ~PNL_RgeHandler();
virtual Bool_t IsValid() { return fIsValid; }
virtual Int_t GetRgeEnergyIndex(const Double_t energy);
virtual Double_t GetRgeValue(const Double_t energy, const Double_t dist);
virtual Double_t GetRgeValue(const Int_t index, const Double_t dist);
private:
Bool_t fIsValid;
PNL_RgeDataList fRgeDataList;
virtual Bool_t LoadRgeData(const PStringVector &rgeDataPathList, const PDoubleVector &rgeDataEnergyList);
};
#endif // _PNL_RGEHANDLER_H_

70
src/external/Nonlocal/PNL_StartupHandler.cpp vendored
View File

@ -32,6 +32,8 @@
#include <iostream> #include <iostream>
#include <fstream> #include <fstream>
#include <boost/filesystem.hpp>
#include "PNL_StartupHandler.h" #include "PNL_StartupHandler.h"
ClassImpQ(PNL_StartupHandler) ClassImpQ(PNL_StartupHandler)
@ -44,46 +46,26 @@ ClassImpQ(PNL_StartupHandler)
*/ */
PNL_StartupHandler::PNL_StartupHandler() PNL_StartupHandler::PNL_StartupHandler()
{ {
fIsValid = true;
fStartupFileFound = false;
fStartupFilePath = "";
fFourierPoints = 0;
fTrimSpDataPath = TString("");
// get default path (for the moment only linux like) // get default path (for the moment only linux like)
char startup_path_name[512]; char startup_path_name[512];
char *home_str=0; char *home_str=0;
// check if the startup file is found in the current directory // check if the startup file is found in the current directory
strcpy(startup_path_name, "./nonlocal_startup.xml"); strcpy(startup_path_name, "./nonlocal_startup.xml");
if (StartupFileExists(startup_path_name)) { if (boost::filesystem::exists(startup_path_name)) {
fStartupFileFound = true; fStartupFileFound = true;
fStartupFilePath = TString(startup_path_name); fStartupFilePath = TString(startup_path_name);
} else { // startup file is not found in the current directory } else { // startup file is not found in the current directory
std::cout << std::endl << "PNL_StartupHandler(): **WARNING** Couldn't find nonlocal_startup.xml in the current directory, will try default one." << std::endl; std::cout << std::endl << "PNL_StartupHandler(): **WARNING** Couldn't find nonlocal_startup.xml in the current directory, will try default one." << std::endl;
home_str = getenv("HOME"); home_str = getenv("HOME");
snprintf(startup_path_name, sizeof(startup_path_name), "%s/.musrfit/nonlocal_startup.xml", home_str); snprintf(startup_path_name, sizeof(startup_path_name), "%s/.musrfit/nonlocal_startup.xml", home_str);
if (StartupFileExists(startup_path_name)) { if (boost::filesystem::exists(startup_path_name)) {
fStartupFileFound = true; fStartupFileFound = true;
fStartupFilePath = TString(startup_path_name); fStartupFilePath = TString(startup_path_name);
} }
} }
} }
//--------------------------------------------------------------------------
// Destructor
//--------------------------------------------------------------------------
/**
*
*/
PNL_StartupHandler::~PNL_StartupHandler()
{
fTrimSpDataPathList.clear();
fTrimSpDataEnergyList.clear();
}
//-------------------------------------------------------------------------- //--------------------------------------------------------------------------
// OnStartDocument // OnStartDocument
//-------------------------------------------------------------------------- //--------------------------------------------------------------------------
@ -124,10 +106,6 @@ void PNL_StartupHandler::OnStartElement(const char *str, const TList *attributes
{ {
if (!strcmp(str, "fourier_points")) { if (!strcmp(str, "fourier_points")) {
fKey = eFourierPoints; fKey = eFourierPoints;
} else if (!strcmp(str, "data_path")) {
fKey = eDataPath;
} else if (!strcmp(str, "energy")) {
fKey = eEnergy;
} }
} }
@ -167,23 +145,6 @@ void PNL_StartupHandler::OnCharacters(const char *str)
std::cout << std::endl; std::cout << std::endl;
} }
break; break;
case eDataPath:
fTrimSpDataPath = str;
break;
case eEnergy:
tstr = str;
if (tstr.IsFloat()) {
fTrimSpDataEnergyList.push_back(tstr.Atof());
tstr = fTrimSpDataPath;
tstr += str;
tstr += ".rge";
fTrimSpDataPathList.push_back(tstr);
} else {
std::cout << std::endl << "PNL_StartupHandler::OnCharacters: **ERROR** when finding energy:";
std::cout << std::endl << "\"" << str << "\" is not a floating point number, will ignore it and use the default value.";
std::cout << std::endl;
}
break;
default: default:
break; break;
} }
@ -257,29 +218,6 @@ void PNL_StartupHandler::OnCdataBlock(const char *str, Int_t len)
// nothing to be done for now // nothing to be done for now
} }
//--------------------------------------------------------------------------
// StartupFileExists
//--------------------------------------------------------------------------
/**
* <p>
*
*/
bool PNL_StartupHandler::StartupFileExists(char *fln)
{
bool result = false;
std::ifstream ifile(fln);
if (ifile.fail()) {
result = false;
} else {
result = true;
ifile.close();
}
return result;
}
// ------------------------------------------------------------------------- // -------------------------------------------------------------------------
// end // end
// ------------------------------------------------------------------------- // -------------------------------------------------------------------------

24
src/external/Nonlocal/PNL_StartupHandler.h vendored
View File

@ -5,8 +5,6 @@
Author: Andreas Suter Author: Andreas Suter
e-mail: andreas.suter@psi.ch e-mail: andreas.suter@psi.ch
$Id$
***************************************************************************/ ***************************************************************************/
/*************************************************************************** /***************************************************************************
@ -36,13 +34,11 @@
#include <TQObject.h> #include <TQObject.h>
#include <TString.h> #include <TString.h>
#include "PNonlocal.h"
class PNL_StartupHandler : public TObject class PNL_StartupHandler : public TObject
{ {
public: public:
PNL_StartupHandler(); PNL_StartupHandler();
virtual ~PNL_StartupHandler(); virtual ~PNL_StartupHandler() {}
virtual void OnStartDocument(); // SLOT virtual void OnStartDocument(); // SLOT
virtual void OnEndDocument(); // SLOT virtual void OnEndDocument(); // SLOT
@ -58,26 +54,18 @@ class PNL_StartupHandler : public TObject
virtual bool IsValid() { return fIsValid; } virtual bool IsValid() { return fIsValid; }
virtual TString GetStartupFilePath() { return fStartupFilePath; } virtual TString GetStartupFilePath() { return fStartupFilePath; }
virtual const Int_t GetFourierPoints() const { return fFourierPoints; } virtual const Int_t GetFourierPoints() const { return fFourierPoints; }
virtual const PStringVector GetTrimSpDataPathList() const { return fTrimSpDataPathList; }
virtual const PDoubleVector GetTrimSpDataVectorList() const { return fTrimSpDataEnergyList; }
virtual bool StartupFileFound() { return fStartupFileFound; } virtual bool StartupFileFound() { return fStartupFileFound; }
private: private:
enum EKeyWords {eEmpty, eComment, eFourierPoints, eDataPath, eEnergy}; enum EKeyWords {eEmpty, eFourierPoints};
EKeyWords fKey; EKeyWords fKey;
bool fIsValid; bool fIsValid{true};
bool fStartupFileFound; bool fStartupFileFound{false};
TString fStartupFilePath; TString fStartupFilePath{""};
Int_t fFourierPoints; Int_t fFourierPoints{0};
TString fTrimSpDataPath;
PStringVector fTrimSpDataPathList;
PDoubleVector fTrimSpDataEnergyList;
bool StartupFileExists(char *fln);
ClassDef(PNL_StartupHandler, 1) ClassDef(PNL_StartupHandler, 1)
}; };

65
src/external/Nonlocal/PNonlocal.h vendored
View File

@ -1,65 +0,0 @@
/***************************************************************************
PNonlocal.h
Author: Andreas Suter
e-mail: andreas.suter@psi.ch
***************************************************************************/
/***************************************************************************
* Copyright (C) 2009-2021 by Andreas Suter *
* andreas.suter@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. *
***************************************************************************/
#ifndef _PNONLOCAL_H_
#define _PNONLOCAL_H_
#include <vector>
#include <TString.h>
//-------------------------------------------------------------
/**
* <p> typedef to make to code more readable.
*/
typedef std::vector<TString> PStringVector;
//-------------------------------------------------------------
/**
* <p> typedef to make to code more readable.
*/
typedef std::vector<Double_t> PDoubleVector;
//-------------------------------------------------------------
/**
* <p>
*/
typedef struct {
Double_t energy;
PDoubleVector stoppingDistance;
PDoubleVector stoppingAmplitude;
} PNL_RgeData;
//-------------------------------------------------------------
/**
* <p>
*/
typedef std::vector<PNL_RgeData> PNL_RgeDataList;
#endif // _PNONLOCAL_H_

32
src/external/Nonlocal/nonlocal_startup.xml vendored
View File

@ -1,25 +1,27 @@
<?xml version="1.0" encoding="UTF-8"?> <?xml version="1.0" encoding="UTF-8"?>
<nonlocal xmlns="http://nemu.web.psi.ch/musrfit/nonlocal"> <nonlocal xmlns="http://nemu.web.psi.ch/musrfit/nonlocal">
<comment> <comment>
trim.sp meta information Fourier and TrimSp information
</comment> </comment>
<nonlocal_par> <nonlocal_par>
<fourier_points>262144</fourier_points> <fourier_points>262144</fourier_points>
</nonlocal_par> </nonlocal_par>
<trim_sp_part> <trim_sp>
<data_path>/afs/psi.ch/project/nemu/analysis/2009/Nonlocal/trimsp/InSne</data_path> <data_path>./profiles/</data_path>
<rge_fln_pre>Sn_E</rge_fln_pre>
<energy_list> <energy_list>
<energy>2.5</energy> <energy>1000</energy>
<energy>4.0</energy> <energy>2000</energy>
<energy>6.0</energy> <energy>4000</energy>
<energy>8.0</energy> <energy>6000</energy>
<energy>10.0</energy> <energy>8000</energy>
<energy>12.5</energy> <energy>10000</energy>
<energy>14.1</energy> <energy>12000</energy>
<energy>17.5</energy> <energy>14100</energy>
<energy>22.0</energy> <energy>18000</energy>
<energy>25.0</energy> <energy>22000</energy>
<energy>28.2</energy> <energy>25000</energy>
<energy>27300</energy>
</energy_list> </energy_list>
</trim_sp_part> </trim_sp>
</nonlocal> </nonlocal>