Changeset 2f9f1ec in sasmodels
- Timestamp:
- Apr 20, 2017 3:57:02 PM (8 years ago)
- Branches:
- costrafo411
- Children:
- 65ceb7d, 7cb1c36
- Parents:
- 2cdc35b (diff), 630156b (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent. - Files:
-
- 13 added
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
example/oriented_usans.py
rea75043 r1cd24b4 19 19 scale=0.08, background=0, 20 20 sld=.291, sld_solvent=7.105, 21 r _polar=1800, r_polar_pd=0.222296, r_polar_pd_n=0,22 r _equatorial=2600, r_equatorial_pd=0.28, r_equatorial_pd_n=0,21 radius_polar=1800, radius_polar_pd=0.222296, radius_polar_pd_n=0, 22 radius_equatorial=2600, radius_equatorial_pd=0.28, radius_equatorial_pd_n=0, 23 23 theta=60, theta_pd=0, theta_pd_n=0, 24 24 phi=60, phi_pd=0, phi_pd_n=0, … … 26 26 27 27 # SET THE FITTING PARAMETERS 28 model.r _polar.range(1000, 10000)29 model.r _equatorial.range(1000, 10000)28 model.radius_polar.range(1000, 10000) 29 model.radius_equatorial.range(1000, 10000) 30 30 model.theta.range(0, 360) 31 31 model.phi.range(0, 360) -
sasmodels/compare.py
r650c6d2 r630156b 73 73 -1d*/-2d computes 1d or 2d data 74 74 -preset*/-random[=seed] preset or random parameters 75 -mono /-poly* force monodisperse/polydisperse75 -mono*/-poly force monodisperse or allow polydisperse demo parameters 76 76 -magnetic/-nonmagnetic* suppress magnetism 77 77 -cutoff=1e-5* cutoff value for including a point in polydispersity … … 84 84 -edit starts the parameter explorer 85 85 -default/-demo* use demo vs default parameters 86 -h tml shows the model docs instead of running the model86 -help/-html shows the model docs instead of running the model 87 87 -title="note" adds note to the plot title, after the model name 88 88 -data="path" uses q, dq from the data file … … 753 753 comp = opts['engines'][1] if have_comp else None 754 754 data = opts['data'] 755 use_data = have_base ^ have_comp755 use_data = (opts['datafile'] is not None) and (have_base ^ have_comp) 756 756 757 757 # Plot if requested 758 758 view = opts['view'] 759 759 import matplotlib.pyplot as plt 760 if limits is None :760 if limits is None and not use_data: 761 761 vmin, vmax = np.Inf, -np.Inf 762 762 if have_base: … … 836 836 'linear', 'log', 'q4', 837 837 'hist', 'nohist', 838 'edit', 'html', 838 'edit', 'html', 'help', 839 839 'demo', 'default', 840 840 ]) … … 947 947 'cutoff' : 0.0, 948 948 'seed' : -1, # default to preset 949 'mono' : False,949 'mono' : True, 950 950 # Default to magnetic a magnetic moment is set on the command line 951 951 'magnetic' : False, … … 958 958 'html' : False, 959 959 'title' : None, 960 'data ': None,960 'datafile' : None, 961 961 } 962 962 engines = [] … … 980 980 elif arg.startswith('-random='): opts['seed'] = int(arg[8:]) 981 981 elif arg.startswith('-title='): opts['title'] = arg[7:] 982 elif arg.startswith('-data='): opts['data '] = arg[6:]982 elif arg.startswith('-data='): opts['datafile'] = arg[6:] 983 983 elif arg == '-random': opts['seed'] = np.random.randint(1000000) 984 984 elif arg == '-preset': opts['seed'] = -1 … … 1005 1005 elif arg == '-default': opts['use_demo'] = False 1006 1006 elif arg == '-html': opts['html'] = True 1007 elif arg == '-help': opts['html'] = True 1007 1008 # pylint: enable=bad-whitespace 1008 1009 … … 1121 1122 1122 1123 # Create the computational engines 1123 if opts['data '] is not None:1124 data = load_data(os.path.expanduser(opts['data ']))1124 if opts['datafile'] is not None: 1125 data = load_data(os.path.expanduser(opts['datafile'])) 1125 1126 else: 1126 1127 data, _ = make_data(opts) -
sasmodels/data.py
ra769b54 r630156b 51 51 from sas.sascalc.dataloader.loader import Loader # type: ignore 52 52 loader = Loader() 53 data = loader.load(filename) 54 if data is None: 53 # Allow for one part in multipart file 54 if '[' in filename: 55 filename, indexstr = filename[:-1].split('[') 56 index = int(indexstr) 57 else: 58 index = None 59 datasets = loader.load(filename) 60 if datasets is None: 55 61 raise IOError("Data %r could not be loaded" % filename) 62 if not isinstance(datasets, list): 63 datasets = [datasets] 64 if index is None and len(datasets) > 1: 65 raise ValueError("Need to specify filename[index] for multipart data") 66 data = datasets[index if index is not None else 0] 56 67 if hasattr(data, 'x'): 57 68 data.qmin, data.qmax = data.x.min(), data.x.max() 58 69 data.mask = (np.isnan(data.y) if data.y is not None 59 70 else np.zeros_like(data.x, dtype='bool')) 71 elif hasattr(data, 'qx_data'): 72 data.mask = ~data.mask 60 73 return data 61 74 … … 450 463 if view is 'log': 451 464 mtheory[mtheory <= 0] = masked 452 plt.plot(data.x, scale*mtheory, '-' , hold=True)465 plt.plot(data.x, scale*mtheory, '-') 453 466 all_positive = all_positive and (mtheory > 0).all() 454 467 some_present = some_present or (mtheory.count() > 0) … … 457 470 plt.ylim(*limits) 458 471 459 plt.xscale('linear' if not some_present or non_positive_x else view) 472 plt.xscale('linear' if not some_present or non_positive_x 473 else view if view is not None 474 else 'log') 460 475 plt.yscale('linear' 461 476 if view == 'q4' or not some_present or not all_positive 462 else view) 477 else view if view is not None 478 else 'log') 463 479 plt.xlabel("$q$/A$^{-1}$") 464 480 plt.ylabel('$I(q)$') 481 title = ("data and model" if use_theory and use_data 482 else "data" if use_data 483 else "model") 484 plt.title(title) 465 485 466 486 if use_calc: … … 482 502 if num_plots > 1: 483 503 plt.subplot(1, num_plots, use_calc + 2) 484 plt.plot(data.x, mresid, ' -')504 plt.plot(data.x, mresid, '.') 485 505 plt.xlabel("$q$/A$^{-1}$") 486 506 plt.ylabel('residuals') 487 plt.xscale('linear' if not some_present or non_positive_x else view) 507 plt.xscale('linear') 508 plt.title('(model - Iq)/dIq') 488 509 489 510 … … 512 533 if theory is not None: 513 534 if is_tof: 514 plt.plot(data.x, np.log(theory)/(data.lam*data.lam), '-' , hold=True)535 plt.plot(data.x, np.log(theory)/(data.lam*data.lam), '-') 515 536 else: 516 plt.plot(data.x, theory, '-' , hold=True)537 plt.plot(data.x, theory, '-') 517 538 if limits is not None: 518 539 plt.ylim(*limits) -
sasmodels/models/onion.py
rc3ccaec rbccb40f 1 1 r""" 2 2 This model provides the form factor, $P(q)$, for a multi-shell sphere where 3 the scattering length density (SLD) of theeach shell is described by an3 the scattering length density (SLD) of each shell is described by an 4 4 exponential, linear, or constant function. The form factor is normalized by 5 5 the volume of the sphere where the SLD is not identical to the SLD of the … … 142 142 143 143 For $A = 0$, the exponential function has no dependence on the radius (so that 144 $\rho_\text{out}$ is ignored this case) and becomes flat. We set the constant144 $\rho_\text{out}$ is ignored in this case) and becomes flat. We set the constant 145 145 to $\rho_\text{in}$ for convenience, and thus the form factor contributed by 146 146 the shells is … … 346 346 # flat shell 347 347 z.append(z_current + thickness[k]) 348 rho.append(sld_ out[k])348 rho.append(sld_in[k]) 349 349 else: 350 350 # exponential shell … … 357 357 z.append(z_current+z_shell) 358 358 rho.append(slope*exp(A[k]*z_shell/thickness[k]) + const) 359 359 360 360 # add in the solvent 361 361 z.append(z[-1]) -
sasmodels/models/spherical_sld.py
r3330bb4 r63a7fe8 199 199 category = "shape:sphere" 200 200 201 SHAPES = [ ["erf(|nu|*z)", "Rpow(z^|nu|)", "Lpow(z^|nu|)",202 "Rexp(-|nu|z)", "Lexp(-|nu|z)"]]201 SHAPES = ["erf(|nu|*z)", "Rpow(z^|nu|)", "Lpow(z^|nu|)", 202 "Rexp(-|nu|z)", "Lexp(-|nu|z)"] 203 203 204 204 # pylint: disable=bad-whitespace, line-too-long … … 209 209 ["thickness[n_shells]", "Ang", 100.0, [0, inf], "volume", "thickness shell"], 210 210 ["interface[n_shells]", "Ang", 50.0, [0, inf], "volume", "thickness of the interface"], 211 ["shape[n_shells]", "", 0, SHAPES,"", "interface shape"],211 ["shape[n_shells]", "", 0, [SHAPES], "", "interface shape"], 212 212 ["nu[n_shells]", "", 2.5, [0, inf], "", "interface shape exponent"], 213 213 ["n_steps", "", 35, [0, inf], "", "number of steps in each interface (must be an odd integer)"], -
sasmodels/rst2html.py
rc4e3215 rf2f5413 155 155 return frame 156 156 157 def view_html_wxapp(html, url=""): 158 import wx # type: ignore 159 app = wx.App() 160 frame = wxview(html, url) 161 app.MainLoop() 162 163 def view_url_wxapp(url): 164 import wx # type: ignore 165 from wx.html2 import WebView 166 app = wx.App() 167 frame = wx.Frame(None, -1, size=(850, 540)) 168 view = WebView.New(frame) 169 view.LoadURL(url) 170 frame.Show() 171 app.MainLoop() 172 157 173 def qtview(html, url=""): 158 174 try: … … 167 183 return helpView 168 184 169 def view_html_wxapp(html, url=""):170 import wx # type: ignore171 app = wx.App()172 frame = wxview(html, url)173 app.MainLoop()174 175 185 def view_html_qtapp(html, url=""): 176 186 import sys … … 183 193 sys.exit(app.exec_()) 184 194 185 def view_rst_app(filename, qt=False): 195 def view_url_qtapp(url): 196 import sys 197 try: 198 from PyQt5.QtWidgets import QApplication 199 except ImportError: 200 from PyQt4.QtGui import QApplication 201 app = QApplication([]) 202 try: 203 from PyQt5.QtWebKitWidgets import QWebView 204 from PyQt5.QtCore import QUrl 205 except ImportError: 206 from PyQt4.QtWebkit import QWebView 207 from PyQt4.QtCore import QUrl 208 frame = QWebView() 209 frame.load(QUrl(url)) 210 frame.show() 211 sys.exit(app.exec_()) 212 213 def view_help(filename, qt=False): 186 214 import os 187 html = load_rst_as_html(filename) 188 url="file://"+os.path.abspath(filename)+"/" 189 if qt: 190 view_html_qtapp(html, url) 215 url="file:///"+os.path.abspath(filename).replace("\\","/") 216 if filename.endswith('.rst'): 217 html = load_rst_as_html(filename) 218 if qt: 219 view_html_qtapp(html, url) 220 else: 221 view_html_wxapp(html, url) 191 222 else: 192 view_html_wxapp(html, url) 223 if qt: 224 view_url_qtapp(url) 225 else: 226 view_url_wxapp(url) 193 227 194 228 if __name__ == "__main__": 195 229 import sys 196 view_ rst_app(sys.argv[1], qt=True)197 230 view_help(sys.argv[1], qt=True) 231 -
sasmodels/direct_model.py
ra769b54 r2cdc35b 202 202 203 203 if self.data_type == 'sesans': 204 q = sesans.make_q(data.sample.zacceptance, data.Rmax) 205 index = slice(None, None) 206 res = None 207 if data.y is not None: 208 Iq, dIq = data.y, data.dy 204 from sas.sascalc.data_util.nxsunit import Converter 205 qmax, qunits = data.sample.zacceptance 206 SElength = Converter(data._xunit)(data.x, "A") 207 zaccept = Converter(qunits)(qmax, "1/A"), 208 Rmax = 10000000 209 index = slice(None, None) # equivalent to index [:] 210 Iq = data.y[index] 211 dIq = data.dy[index] 212 oriented = getattr(data, 'oriented', False) 213 if oriented: 214 res = sesans.OrientedSesansTransform(data.x[index], SElength, zaccept, Rmax) 215 # Oriented sesans transform produces q_calc = [qx, qy] 216 q_vectors = res.q_calc 209 217 else: 210 Iq, dIq = None, None 211 #self._theory = np.zeros_like(q) 212 q_vectors = [q] 213 q_mono = sesans.make_all_q(data) 218 res = sesans.SesansTransform(data.x[index], SElength, zaccept, Rmax) 219 # Unoriented sesans transform produces q_calc = q 220 q_vectors = [res.q_calc] 214 221 elif self.data_type == 'Iqxy': 215 222 #if not model.info.parameters.has_2d: … … 230 237 #self._theory = np.zeros_like(self.Iq) 231 238 q_vectors = res.q_calc 232 q_mono = []233 239 elif self.data_type == 'Iq': 234 240 index = (data.x >= data.qmin) & (data.x <= data.qmax) … … 255 261 #self._theory = np.zeros_like(self.Iq) 256 262 q_vectors = [res.q_calc] 257 q_mono = []258 263 elif self.data_type == 'Iq-oriented': 259 264 index = (data.x >= data.qmin) & (data.x <= data.qmax) … … 272 277 qy_width=data.dxl[index]) 273 278 q_vectors = res.q_calc 274 q_mono = []275 279 else: 276 280 raise ValueError("Unknown data type") # never gets here … … 279 283 # so we can save/restore state 280 284 self._kernel_inputs = q_vectors 281 self._kernel_mono_inputs = q_mono282 285 self._kernel = None 283 286 self.Iq, self.dIq, self.index = Iq, dIq, index … … 306 309 if self._kernel is None: 307 310 self._kernel = self._model.make_kernel(self._kernel_inputs) 308 self._kernel_mono = (309 self._model.make_kernel(self._kernel_mono_inputs)310 if self._kernel_mono_inputs else None)311 311 312 312 Iq_calc = call_kernel(self._kernel, pars, cutoff=cutoff) … … 316 316 # TODO: refactor so we don't store the result in the model 317 317 self.Iq_calc = None 318 if self.data_type == 'sesans': 319 Iq_mono = (call_kernel(self._kernel_mono, pars, mono=True) 320 if self._kernel_mono_inputs else None) 321 result = sesans.transform(self._data, 322 self._kernel_inputs[0], Iq_calc, 323 self._kernel_mono_inputs, Iq_mono) 324 else: 325 result = self.resolution.apply(Iq_calc) 326 if hasattr(self.resolution, 'nx'): 327 self.Iq_calc = ( 328 self.resolution.qx_calc, self.resolution.qy_calc, 329 np.reshape(Iq_calc, (self.resolution.ny, self.resolution.nx)) 330 ) 318 result = self.resolution.apply(Iq_calc) 319 if hasattr(self.resolution, 'nx'): 320 self.Iq_calc = ( 321 self.resolution.qx_calc, self.resolution.qy_calc, 322 np.reshape(Iq_calc, (self.resolution.ny, self.resolution.nx)) 323 ) 331 324 return result 332 325 -
sasmodels/sesans.py
r94d13f1 r2cdc35b 20 20 Spin-Echo SANS transform calculator. Similar to a resolution function, 21 21 the SesansTransform object takes I(q) for the set of *q_calc* values and 22 produces a transformed dataset 22 produces a transformed dataset. 23 23 24 24 *SElength* (A) is the set of spin-echo lengths in the measured data. … … 48 48 49 49 def apply(self, Iq): 50 # tye: (np.ndarray) -> np.ndarray51 50 G0 = np.dot(self._H0, Iq) 52 51 G = np.dot(self._H.T, Iq) … … 73 72 self.q_calc = q 74 73 self._H, self._H0 = H, H0 74 75 class OrientedSesansTransform(object): 76 """ 77 Oriented Spin-Echo SANS transform calculator. Similar to a resolution 78 function, the OrientedSesansTransform object takes I(q) for the set 79 of *q_calc* values and produces a transformed dataset. 80 81 *SElength* (A) is the set of spin-echo lengths in the measured data. 82 83 *zaccept* (1/A) is the maximum acceptance of scattering vector in the spin 84 echo encoding dimension (for ToF: Q of min(R) and max(lam)). 85 86 *Rmax* (A) is the maximum size sensitivity; larger radius requires more 87 computation time. 88 """ 89 #: SElength from the data in the original data units; not used by transform 90 #: but the GUI uses it, so make sure that it is present. 91 q = None # type: np.ndarray 92 93 #: q values to calculate when computing transform 94 q_calc = None # type: np.ndarray 95 96 # transform arrays 97 _cosmat = None # type: np.ndarray 98 _cos0 = None # type: np.ndarray 99 _Iq_shape = None # type: Tuple[int, int] 100 101 def __init__(self, z, SElength, zaccept, Rmax): 102 # type: (np.ndarray, float, float) -> None 103 #import logging; logging.info("creating SESANS transform") 104 self.q = z 105 self._set_cosmat(SElength, zaccept, Rmax) 106 107 def apply(self, Iq): 108 dq = self.q_calc[0][0] 109 Iq = np.reshape(Iq, self._Iq_shape) 110 G0 = self._cos0 * np.sum(Iq) * dq 111 G = np.sum(np.dot(Iq, self._cosmat), axis=0) * dq 112 P = G - G0 113 return P 114 115 def _set_cosmat(self, SElength, zaccept, Rmax): 116 # type: (np.ndarray, float, float) -> None 117 # Force float32 arrays, otherwise run into memory problems on some machines 118 SElength = np.asarray(SElength, dtype='float32') 119 120 # Rmax = #value in text box somewhere in FitPage? 121 q_max = 2 * pi / (SElength[1] - SElength[0]) 122 q_min = 0.1 * 2 * pi / (np.size(SElength) * SElength[-1]) 123 q_min *= 100 124 125 q = np.arange(q_min, q_max, q_min, dtype='float32') 126 dq = q_min 127 128 cos0 = np.float32(dq / (2 * pi)) 129 cosmat = np.float32(dq / (2 * pi)) * np.cos(q[:, None] * SElength[None, :]) 130 131 qx, qy = np.meshgrid(q, q) 132 self._Iq_shape = qx.shape 133 self.q_calc = qx.flatten(), qy.flatten() 134 self._cosmat, self._cos0 = cosmat, cos0
Note: See TracChangeset
for help on using the changeset viewer.