source: sasview/src/sas/sascalc/corfunc/transform_thread.py @ f7e6b30

ESS_GUIESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_opencl
Last change on this file since f7e6b30 was f7e6b30, checked in by awashington, 5 years ago

Remote artificial limits on corfunc range.

With the new navigation toolbars, they aren't needed.

  • Property mode set to 100644
File size: 3.7 KB
RevLine 
[a2db1ab]1from sas.sascalc.data_util.calcthread import CalcThread
2from sas.sascalc.dataloader.data_info import Data1D
3from scipy.fftpack import dct
[5e2267e]4from scipy.integrate import trapz, cumtrapz
[a2db1ab]5import numpy as np
6from time import sleep
7
[d03228e]8class FourierThread(CalcThread):
[c728295]9    def __init__(self, raw_data, extrapolated_data, bg, updatefn=None,
10        completefn=None):
[a2db1ab]11        CalcThread.__init__(self, updatefn=updatefn, completefn=completefn)
12        self.data = raw_data
13        self.background = bg
14        self.extrapolation = extrapolated_data
15
[37d53225]16    def check_if_cancelled(self):
17        if self.isquit():
18            self.update("Fourier transform cancelled.")
19            self.complete(transforms=None)
20            return True
21        return False
22
[a2db1ab]23    def compute(self):
24        qs = self.extrapolation.x
25        iqs = self.extrapolation.y
26        q = self.data.x
27        background = self.background
28
[412c509]29        xs = np.pi*np.arange(len(qs),dtype=np.float32)/(q[1]-q[0])/len(qs)
30
[a2db1ab]31        self.ready(delay=0.0)
[412c509]32        self.update(msg="Fourier transform in progress.")
[a2db1ab]33        self.ready(delay=0.0)
[412c509]34
[37d53225]35        if self.check_if_cancelled(): return
[a2db1ab]36        try:
[37d53225]37            # ----- 1D Correlation Function -----
[457f735]38            gamma1 = dct((iqs-background)*qs**2)
[9b90bf8]39            Q = gamma1.max()
40            gamma1 /= Q
[412c509]41
[37d53225]42            if self.check_if_cancelled(): return
43
44            # ----- 3D Correlation Function -----
[412c509]45            # gamma3(R) = 1/R int_{0}^{R} gamma1(x) dx
[2a54ba5]46            # trapz uses the trapezium rule to calculate the integral
[a859f99]47            # gamma3 = [trapz(gamma1[:n], xs[:n])/xs[n-1] for n in range(2, len(xs[mask]) + 1)]j
48            # gamma3.insert(0, 1.0) # Gamma_3(0) is defined as 1
[f7e6b30]49            n = len(xs)
[a859f99]50            gamma3 = cumtrapz(gamma1[:n], xs[:n])/xs[1:n]
[5e2267e]51            gamma3 = np.hstack((1.0, gamma3)) # Gamma_3(0) is defined as 1
[37d53225]52
53            if self.check_if_cancelled(): return
54
55            # ----- Interface Distribution function -----
[a309667]56            idf = dct(-qs**4 * (iqs-background))
[37d53225]57
[a309667]58            if self.check_if_cancelled(): return
[37d53225]59
[a309667]60            # Manually calculate IDF(0.0), since scipy DCT tends to give us a
61            # very large negative value.
62            # IDF(x) = int_0^inf q^4 * I(q) * cos(q*x) * dq
63            # => IDF(0) = int_0^inf q^4 * I(q) * dq
[9b90bf8]64            idf[0] = trapz(-qs**4 * (iqs-background), qs)
[a309667]65            idf /= Q # Normalise using scattering invariant
[37d53225]66
[412c509]67        except Exception as e:
68            import logging
69            logger = logging.getLogger(__name__)
70            logger.error(e)
71
[a2db1ab]72            self.update(msg="Fourier transform failed.")
[412c509]73            self.complete(transforms=None)
[a2db1ab]74            return
75        if self.isquit():
76            return
77        self.update(msg="Fourier transform completed.")
78
[457f735]79        transform1 = Data1D(xs, gamma1)
[f7e6b30]80        transform3 = Data1D(xs, gamma3)
[9b90bf8]81        idf = Data1D(xs, idf)
[a2db1ab]82
[37d53225]83        transforms = (transform1, transform3, idf)
[457f735]84
85        self.complete(transforms=transforms)
[d03228e]86
87class HilbertThread(CalcThread):
88    def __init__(self, raw_data, extrapolated_data, bg, updatefn=None,
89        completefn=None):
90        CalcThread.__init__(self, updatefn=updatefn, completefn=completefn)
91        self.data = raw_data
92        self.background = bg
93        self.extrapolation = extrapolated_data
94
95    def compute(self):
96        qs = self.extrapolation.x
97        iqs = self.extrapolation.y
98        q = self.data.x
99        background = self.background
100
101        self.ready(delay=0.0)
102        self.update(msg="Starting Hilbert transform.")
103        self.ready(delay=0.0)
104        if self.isquit():
105            return
106
107        # TODO: Implement hilbert transform
108
109        self.update(msg="Hilbert transform completed.")
110
[626f833]111        self.complete(transforms=None)
Note: See TracBrowser for help on using the repository browser.