Compare commits
42 Commits
Author | SHA1 | Date | |
---|---|---|---|
d1c2632e81 | |||
e19795b51c | |||
9c56171524 | |||
6980df1599 | |||
f6b7ebec07 | |||
8d1ccd22fa | |||
0b52fef549 | |||
b6b98d292a | |||
d90959c6b6 | |||
e459bd5e54 | |||
24f20f8850 | |||
fc91bf83fe | |||
86f285fba5 | |||
fa84b0382e | |||
4b75aa9267 | |||
f185b24967 | |||
d07b85ae27 | |||
7ad1e4b843 | |||
ff2ff01da7 | |||
d9f1c0b8c2 | |||
9039c44ce7 | |||
ac6b734f81 | |||
9babb73f3a | |||
4f0a7827ba | |||
05862730a0 | |||
7fe89eff7f | |||
f94f78893c | |||
fda3257424 | |||
f30ff3b758 | |||
7b61c1244d | |||
e0c287d8a9 | |||
a8fcd658d9 | |||
dd471ae294 | |||
ab586ac39a | |||
b355aab99d | |||
03cdc225ca | |||
e87c6bf2c1 | |||
cc7572fe14 | |||
ef66cf584a | |||
f0448fac0f | |||
749a78b550 | |||
50a811b7ec |
@ -133,7 +133,7 @@ class QAsciiReader(QtWidgets.QDialog, Ui_ascii_reader):
|
||||
if self.column_checkBox.isChecked() and self.line_spinBox.isEnabled():
|
||||
header_line = self.reader.header[self.line_spinBox.value()-1]
|
||||
header_line = header_line.strip('\n\t\r, ')
|
||||
header_line = re.sub(r'[\t, ;]+(?!\w*})', ';', header_line)
|
||||
header_line = re.sub(r'[\t ;,]+', ';', header_line)
|
||||
|
||||
self.ascii_table.setHorizontalHeaderLabels(header_line.split(';'))
|
||||
|
||||
|
@ -1,16 +0,0 @@
|
||||
/* integrands used in quadrature integration with scipy's LowLevelCallables */
|
||||
#include <math.h>
|
||||
|
||||
double anistropicDiffusion(double x, void *user_data) {
|
||||
double *c = (double *)user_data;
|
||||
|
||||
double q = c[0];
|
||||
double t = c[1];
|
||||
double d_perp = c[2];
|
||||
double d_par = c[3];
|
||||
|
||||
double cos_theta = cos(x);
|
||||
double sin_theta = sin(x);
|
||||
|
||||
return exp(-q * q * t * (d_par * cos_theta * cos_theta + d_perp * sin_theta * sin_theta)) * sin_theta;
|
||||
}
|
Binary file not shown.
@ -5,17 +5,6 @@ from ctypes import CDLL, c_double, c_void_p
|
||||
from ..lib.logger import logger
|
||||
|
||||
|
||||
diffusion_lib = None
|
||||
try:
|
||||
diffusion_lib = CDLL(str(Path(__file__).parents[1] / 'clib' / 'diffusion.so'))
|
||||
|
||||
diffusion_lib.anistropicDiffusion.restype = c_double
|
||||
diffusion_lib.anistropicDiffusion.argtypes = (c_double, c_void_p)
|
||||
|
||||
HAS_C_FUNCS = True
|
||||
except OSError:
|
||||
HAS_C_FUNCS = False
|
||||
|
||||
lib = None
|
||||
try:
|
||||
lib = CDLL(str(Path(__file__).parents[1] / 'clib' / 'integrate.so'))
|
||||
@ -50,8 +39,10 @@ try:
|
||||
lib.energyDistSuscImag.restype = c_double
|
||||
lib.energyDistSuscImag.argtypes = (c_double, c_void_p)
|
||||
|
||||
|
||||
HAS_C_FUNCS = True
|
||||
logger.info('Use C functions')
|
||||
except OSError:
|
||||
HAS_C_FUNCS = False
|
||||
logger.info('Use python functions')
|
||||
|
||||
|
@ -49,7 +49,7 @@ class AsciiReader:
|
||||
with self.fname.open('r') as f:
|
||||
for i, line in enumerate(islice(f, len(self.header)+len(self.lines), num_lines)):
|
||||
line = line.strip('\n\t\r, ')
|
||||
line = re.sub(r'[\t, ;]+(?!\w*})', ';', line)
|
||||
line = re.sub(r'[\t ;,]+', ';', line)
|
||||
line = line.split(';')
|
||||
|
||||
try:
|
||||
@ -146,11 +146,10 @@ class AsciiReader:
|
||||
raw_data = raw_data.reshape((1, *raw_data.shape))
|
||||
|
||||
if len(x) == 0 or raw_data.shape[2] == 1:
|
||||
raw_data = raw_data.reshape(raw_data.shape[0], raw_data.shape[2], raw_data.shape[1])
|
||||
# _temp = np.zeros((raw_data.shape[0], raw_data.shape[2], raw_data.shape[1]))
|
||||
# _temp[:, :, 0] = np.arange(raw_data.shape[1])
|
||||
# _temp[:, :, 1:] = raw_data
|
||||
# raw_data = _temp
|
||||
_temp = np.zeros((raw_data.shape[0], raw_data.shape[1], raw_data.shape[2]+1))
|
||||
_temp[:, :, 0] = np.arange(raw_data.shape[1])
|
||||
_temp[:, :, 1:] = raw_data
|
||||
raw_data = _temp
|
||||
|
||||
if y:
|
||||
y = [i+1 for i in y]
|
||||
|
@ -302,26 +302,15 @@ class HdfReader(HdfNode):
|
||||
|
||||
def make_signal(self, node, flag: str = 'fid', value: str = None, group: str = None):
|
||||
if value is None:
|
||||
data_name = node.name
|
||||
value = self._get_parameter_values(node, node.parameter)
|
||||
else:
|
||||
try:
|
||||
data_name = f"{value}={node.parameter[value]}"
|
||||
value = node.parameter[value]
|
||||
except KeyError:
|
||||
print(node.title_parameter)
|
||||
try:
|
||||
temp = node
|
||||
while value != temp.title_parameter[0][0]:
|
||||
if temp.parent is None:
|
||||
break
|
||||
temp = temp.parent
|
||||
|
||||
value = temp.title_parameter[0][1]
|
||||
data_name = temp.name
|
||||
value = node.title_parameter[1][value]
|
||||
except KeyError:
|
||||
print(f'{value} is not a valid key for {node.name}')
|
||||
data_name = node.name
|
||||
value = None
|
||||
|
||||
if group is None:
|
||||
@ -354,11 +343,11 @@ class HdfReader(HdfNode):
|
||||
dw = float(index['dwelltime'])
|
||||
if flag == 'fid':
|
||||
x = np.arange(len(y)) * dw
|
||||
ret = FID(x, y, name=data_name, value=value, group=group, filename=self.file.filename)
|
||||
ret = FID(x, y, name=node.name, value=value, group=group, filename=self.file.filename)
|
||||
|
||||
elif flag == 'spectrum':
|
||||
x = np.linspace(-1/dw, 1/dw, num=len(y))
|
||||
ret = Spectrum(x, y, name=data_name, value=value, group=group, filename=self.file.filename)
|
||||
ret = Spectrum(x, y, name=node.name, value=value, group=group, filename=self.file.filename)
|
||||
else:
|
||||
raise ValueError(f'{flag} unknown, use `fid` or `spectrum`.')
|
||||
|
||||
|
@ -203,31 +203,6 @@ class Sinc:
|
||||
return c * np.sinc(((x-x0)/w)/np.pi)
|
||||
|
||||
|
||||
class Sigmoid:
|
||||
type = 'Basic'
|
||||
name = 'Sigmoid'
|
||||
equation = 'C / [1 + exp(-a * (x - x_{0})] + y_{0}'
|
||||
params = ['C', 'a', 'x_{0}', 'y_{0}']
|
||||
|
||||
@staticmethod
|
||||
def func(x, c, a, x0, y0):
|
||||
"""
|
||||
Sigmoid function
|
||||
|
||||
.. math::
|
||||
y = C / [1 + exp(-a * (x - x_0))] + y_0
|
||||
|
||||
Args:
|
||||
x (array_like): Input values
|
||||
c (float): Prefactor
|
||||
a (float): Steepness of the sigmoid
|
||||
x0 (float): x position of the sigmoid's midpoint
|
||||
y0 (float): y position of the sigmoid's midpoint
|
||||
|
||||
"""
|
||||
return c / (1 + np.exp(-a * (x - x0))) + y0
|
||||
|
||||
|
||||
class Sine:
|
||||
"""
|
||||
Wavy sine function
|
||||
|
@ -1,11 +1,7 @@
|
||||
from ctypes import c_double, cast, c_void_p, pointer
|
||||
|
||||
import numpy as np
|
||||
from scipy import special as special, LowLevelCallable
|
||||
from scipy.integrate import quad
|
||||
from scipy import special as special
|
||||
|
||||
from ..utils import gamma
|
||||
from nmreval.distributions.helper import HAS_C_FUNCS, diffusion_lib
|
||||
|
||||
|
||||
class Diffusion:
|
||||
@ -107,36 +103,20 @@ class AnisotropicDiffusion(object):
|
||||
tp = x
|
||||
relax = np.exp(-(tp/trel)**brel)*np.exp(-(tp/trel)**brel)
|
||||
|
||||
q = g * nucleus * tp
|
||||
q_squared = np.power(g * nucleus * tp, 2)
|
||||
t = 2 * tp / 3 + tm
|
||||
z = np.sqrt(q_squared * (d_par - d_perp) * t)
|
||||
# Callaghan eq (6.89)
|
||||
if HAS_C_FUNCS:
|
||||
# divide by 2 to normalize by integral sin(x), x=0..pi
|
||||
diffusion_decay = AnisotropicDiffusion._integrate_c(q, t, d_perp, d_par) / 2
|
||||
else:
|
||||
z = np.sqrt(q**2 * (d_par - d_perp) * t)
|
||||
diffusion_decay = np.exp(-q**2 * t * d_perp) * special.erf(z) / z
|
||||
diffs = np.exp(-q_squared*t*d_perp) * special.erf(z) / z
|
||||
|
||||
return m0 * diffusion_decay * relax
|
||||
|
||||
@staticmethod
|
||||
def _integrate_c(q, t, d_perp, d_par) -> np.ndarray:
|
||||
diffusion_decay = np.zeros_like(t)
|
||||
|
||||
for (i, t_i) in enumerate(t):
|
||||
c = (c_double * 4)(q, t_i, d_perp, d_par)
|
||||
user_data = cast(pointer(c), c_void_p)
|
||||
|
||||
diffusion_decay[i] = quad(LowLevelCallable(diffusion_lib.anistropicDiffusion, user_data), 0, np.pi, epsabs=1e-13)[0]
|
||||
|
||||
return diffusion_decay
|
||||
return m0 * diffs * relax
|
||||
|
||||
|
||||
class Peschier:
|
||||
name = 'Diffusion + Cross-Relaxation'
|
||||
type = 'Diffusion'
|
||||
equation = r'Diffusion with cross-relax f(ast) \rightarrow s(low)'
|
||||
params = ['M_{0}', 'D', 'T_{1f}', 'T_{1s}', 'k', 'p_{f}', 't_{ev}', 'g']
|
||||
params = ['M_{0}', 'D', 'T_{1,f}', 'T_{1,s}', 'k', 'p_{f}', 't_{ev}', 'g']
|
||||
bounds = [(0, None), (0, None), (0, None), (0, None), (0, None), (0, None)]
|
||||
choices = [(r'\gamma', 'nucleus', gamma)]
|
||||
|
||||
|
@ -75,7 +75,7 @@ class TwoSatRecAbsolute:
|
||||
type = 'Relaxation'
|
||||
name = 'Two-step relaxation (abs. int)'
|
||||
equation = r'M_{0} + \Sigma \DeltaM_{i}(1-exp(-(x/T_{1,i})^{\beta_{i}}))'
|
||||
params = [r'\DeltaM_{1}', 'T_{11}', r'\beta_{1}', r'\DeltaM_{2}', 'T_{12}', r'\beta_{2}', 'M_{0}']
|
||||
params = [r'\DeltaM_{1}', 'T_{1,1}', r'\beta_{1}', r'\DeltaM_{2}', 'T_{1,2}', r'\beta_{2}', 'M_{0}']
|
||||
choices = [('Type', 'is_inv', {'Saturation': False, 'Inversion': True})]
|
||||
bounds = [(0, None), (0, None), (0, 1), (0, None), (0, None), (0, 1), (None, None)]
|
||||
|
||||
@ -92,7 +92,7 @@ class TwoSatRecRelative:
|
||||
name = 'Two-step relaxation (rel. int)'
|
||||
equation = r'M_{0} + \DeltaM[R(1-exp(-(x/T_{1,1})^{\beta_{1}})) + \n'\
|
||||
r'(1-R)(1-exp(-(x/T_{1,2})^{\beta_{2}}))]'
|
||||
params = [r'\DeltaM', 'M_{0}', 'T_{11}', r'\beta_{1}', 'T_{12}', r'\beta_{2}', 'R']
|
||||
params = [r'\DeltaM', 'M_{0}', 'T_{1,1}', r'\beta_{1}', 'T_{1,2}', r'\beta_{2}', 'R']
|
||||
choices = [('Type', 'kind', {'Saturation': 'sat', 'Inversion': 'inv'})]
|
||||
bounds = [(0, None), (None, None), (0, None), (0, 1), (0, None), (0, 1), (0, 1)]
|
||||
|
||||
|
Reference in New Issue
Block a user