1 | """ |
---|
2 | Conversion of scattering cross section from SANS (I(q), or rather, ds/dO) in absolute |
---|
3 | units (cm-1)into SESANS correlation function G using a Hankel transformation, then converting |
---|
4 | the SESANS correlation function into polarisation from the SESANS experiment |
---|
5 | |
---|
6 | Everything is in units of metres except specified otherwise (NOT TRUE!!!) |
---|
7 | Everything is in conventional units (nm for spin echo length) |
---|
8 | |
---|
9 | Wim Bouwman (w.g.bouwman@tudelft.nl), June 2013 |
---|
10 | """ |
---|
11 | |
---|
12 | from __future__ import division |
---|
13 | |
---|
14 | import numpy as np # type: ignore |
---|
15 | from numpy import pi, exp # type: ignore |
---|
16 | from scipy.special import jv as besselj |
---|
17 | from scipy.special import j0 |
---|
18 | #from mpmath import j0 |
---|
19 | from mpmath import besselj |
---|
20 | from mpmath import mpf |
---|
21 | from src.sas.sascalc.data_util.nxsunit import Converter |
---|
22 | #from sas.sasgui.perspectives.fitting.fitpage import FitPage |
---|
23 | #import direct_model.DataMixin as model |
---|
24 | |
---|
25 | def make_q(q_max, Rmax): |
---|
26 | r""" |
---|
27 | Return a $q$ vector suitable for SESANS covering from $2\pi/ (10 R_{\max})$ |
---|
28 | to $q_max$. This is the integration range of the Hankel transform; bigger range and |
---|
29 | more points makes a better numerical integration. |
---|
30 | Smaller q_min will increase reliable spin echo length range. |
---|
31 | Rmax is the "radius" of the largest expected object and can be set elsewhere. |
---|
32 | q_max is determined by the acceptance angle of the SESANS instrument. |
---|
33 | """ |
---|
34 | from sas.sascalc.data_util.nxsunit import Converter |
---|
35 | |
---|
36 | q_min = dq = 0.1 * 2*pi / Rmax |
---|
37 | return np.arange(q_min, |
---|
38 | Converter(q_max[1])(q_max[0], |
---|
39 | units="1/A"), |
---|
40 | dq) |
---|
41 | |
---|
42 | def Hankelconstructor(data): |
---|
43 | Rmax = 1000000 |
---|
44 | q_calc = make_q(data.sample.zacceptance, Rmax) |
---|
45 | SElength = Converter(data._xunit)(data.x, "A") |
---|
46 | dq = q_calc[1] - q_calc[0] |
---|
47 | H0 = dq / (2 * pi) * q_calc |
---|
48 | repSE, repq = np.meshgrid(SElength,q_calc) |
---|
49 | H = dq / (2 * pi) * j0(repSE*repq)*repq |
---|
50 | return H0, H, q_calc |
---|
51 | |
---|
52 | def hankeltrafo(H0, H, Iq_calc): |
---|
53 | G0 = np.dot(H0, Iq_calc) |
---|
54 | G = np.dot(H.T, Iq_calc) |
---|
55 | P = G - G0 |
---|
56 | return P # This is the logarithmic Polarization, not the linear one as found in Andersson2008! |
---|