Changeset ed82794 in sasmodels
- Timestamp:
- Jan 29, 2016 4:03:21 AM (9 years ago)
- Branches:
- master, core_shell_microgels, costrafo411, magnetic_model, release_v0.94, release_v0.95, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
- Children:
- 3a45c2c, 504abee
- Parents:
- 168052c (diff), 190fc2b (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:
-
- 24 edited
Legend:
- Unmodified
- Added
- Removed
-
extra/pylint.rc
r37a7252 rd4666ca 7 7 # pygtk.require(). 8 8 #init-hook='import os, sys; sys.path.extend([os.getcwd(),os.path.join(os.getcwd(),"extra")])' 9 init-hook='import sys; sys.path.append('extra') 9 10 10 11 # Profiled execution. … … 64 65 star-args, 65 66 unbalanced-tuple-unpacking, 67 locally-enabled, 66 68 locally-disabled, 67 69 … … 107 109 good-names=_, 108 110 input, 109 i,j,k,n,x,y,z,111 h,i,j,k,n,w,x,y,z, 110 112 q,qx,qy,qz, 111 113 dt,dx,dy,dz,id, 112 Iq,dIq,Qx,Qy,Qz, 114 Iq,dIq,Qx,Qy,Qz,Iqxy, 113 115 p, 114 116 ER, call_ER, VR, call_VR, 117 Rmax, SElength, 115 118 116 119 # Bad variable names which should always be refused, separated by a comma … … 335 338 336 339 # Maximum number of arguments for function / method 337 max-args= 5340 max-args=15 338 341 339 342 # Argument names that match this expression will be ignored. Default to name … … 342 345 343 346 # Maximum number of locals for function / method body 344 max-locals= 15347 max-locals=25 345 348 346 349 # Maximum number of return / yield for function / method body -
sasmodels/bumps_model.py
r37a7252 r190fc2b 12 12 """ 13 13 14 import warnings 15 16 import numpy as np 17 18 from .data import plot_theory 19 from .direct_model import DataMixin 20 14 21 __all__ = [ 15 22 "Model", "Experiment", 16 23 ] 17 18 import warnings19 20 import numpy as np21 22 from .data import plot_theory23 from .direct_model import DataMixin24 24 25 25 # CRUFT: old style bumps wrapper which doesn't separate data and model -
sasmodels/compare.py
rd15a908 r190fc2b 28 28 29 29 from __future__ import print_function 30 31 import sys 32 import math 33 from os.path import basename, dirname, join as joinpath 34 import glob 35 import datetime 36 import traceback 37 38 import numpy as np 39 40 from . import core 41 from . import kerneldll 42 from . import generate 43 from .data import plot_theory, empty_data1D, empty_data2D 44 from .direct_model import DirectModel 45 from .convert import revert_model, constrain_new_to_old 30 46 31 47 USAGE = """ … … 80 96 + USAGE) 81 97 82 83 84 import sys 85 import math 86 from os.path import basename, dirname, join as joinpath 87 import glob 88 import datetime 89 import traceback 90 91 import numpy as np 92 98 kerneldll.ALLOW_SINGLE_PRECISION_DLLS = True 99 100 # List of available models 93 101 ROOT = dirname(__file__) 94 sys.path.insert(0, ROOT) # Make sure sasmodels is first on the path95 96 97 from . import core98 from . import kerneldll99 from . import generate100 from .data import plot_theory, empty_data1D, empty_data2D101 from .direct_model import DirectModel102 from .convert import revert_model, constrain_new_to_old103 kerneldll.ALLOW_SINGLE_PRECISION_DLLS = True104 105 # List of available models106 102 MODELS = [basename(f)[:-3] 107 103 for f in sorted(glob.glob(joinpath(ROOT, "models", "[a-zA-Z]*.py")))] … … 117 113 return dt.total_seconds() 118 114 115 116 class push_seed(object): 117 """ 118 Set the seed value for the random number generator. 119 120 When used in a with statement, the random number generator state is 121 restored after the with statement is complete. 122 123 :Parameters: 124 125 *seed* : int or array_like, optional 126 Seed for RandomState 127 128 :Example: 129 130 Seed can be used directly to set the seed:: 131 132 >>> from numpy.random import randint 133 >>> push_seed(24) 134 <...push_seed object at...> 135 >>> print(randint(0,1000000,3)) 136 [242082 899 211136] 137 138 Seed can also be used in a with statement, which sets the random 139 number generator state for the enclosed computations and restores 140 it to the previous state on completion:: 141 142 >>> with push_seed(24): 143 ... print(randint(0,1000000,3)) 144 [242082 899 211136] 145 146 Using nested contexts, we can demonstrate that state is indeed 147 restored after the block completes:: 148 149 >>> with push_seed(24): 150 ... print(randint(0,1000000)) 151 ... with push_seed(24): 152 ... print(randint(0,1000000,3)) 153 ... print(randint(0,1000000)) 154 242082 155 [242082 899 211136] 156 899 157 158 The restore step is protected against exceptions in the block:: 159 160 >>> with push_seed(24): 161 ... print(randint(0,1000000)) 162 ... try: 163 ... with push_seed(24): 164 ... print(randint(0,1000000,3)) 165 ... raise Exception() 166 ... except: 167 ... print("Exception raised") 168 ... print(randint(0,1000000)) 169 242082 170 [242082 899 211136] 171 Exception raised 172 899 173 """ 174 def __init__(self, seed=None): 175 self._state = np.random.get_state() 176 np.random.seed(seed) 177 178 def __enter__(self): 179 return None 180 181 def __exit__(self, *args): 182 np.random.set_state(self._state) 119 183 120 184 def tic(): … … 179 243 return [0, (2*v if v > 0 else 1)] 180 244 245 181 246 def _randomize_one(p, v): 182 247 """ … … 188 253 return np.random.uniform(*parameter_range(p, v)) 189 254 255 190 256 def randomize_pars(pars, seed=None): 191 257 """ … … 197 263 :func:`constrain_pars` needs to be called afterward.. 198 264 """ 199 np.random.seed(seed)200 # Note: the sort guarantees order `of calls to random number generator201 pars = dict((p, _randomize_one(p, v))202 for p, v in sorted(pars.items()))265 with push_seed(seed): 266 # Note: the sort guarantees order `of calls to random number generator 267 pars = dict((p, _randomize_one(p, v)) 268 for p, v in sorted(pars.items())) 203 269 return pars 204 270 … … 463 529 if Nbase > 0: 464 530 if Ncomp > 0: plt.subplot(131) 465 plot_theory(data, base_value, view=view, plot_data=False, limits=limits)531 plot_theory(data, base_value, view=view, use_data=False, limits=limits) 466 532 plt.title("%s t=%.1f ms"%(base.engine, base_time)) 467 533 #cbar_title = "log I" 468 534 if Ncomp > 0: 469 535 if Nbase > 0: plt.subplot(132) 470 plot_theory(data, comp_value, view=view, plot_data=False, limits=limits)536 plot_theory(data, comp_value, view=view, use_data=False, limits=limits) 471 537 plt.title("%s t=%.1f ms"%(comp.engine, comp_time)) 472 538 #cbar_title = "log I" … … 478 544 err, errstr, errview = abs(relerr), "rel err", "log" 479 545 #err,errstr = base/comp,"ratio" 480 plot_theory(data, None, resid=err, view=errview, plot_data=False)546 plot_theory(data, None, resid=err, view=errview, use_data=False) 481 547 plt.title("max %s = %.3g"%(errstr, max(abs(err)))) 482 548 #cbar_title = errstr if errview=="linear" else "log "+errstr -
sasmodels/compare_many.py
rd15a908 r4f2478e 261 261 262 262 if __name__ == "__main__": 263 #from .compare import push_seed 264 #with push_seed(1): main() 263 265 main() -
sasmodels/core.py
rd15a908 r190fc2b 2 2 Core model handling routines. 3 3 """ 4 __all__ = [5 "list_models", "load_model_definition", "precompile_dll",6 "load_model", "make_kernel", "call_kernel", "call_ER", "call_VR",7 ]8 4 9 5 from os.path import basename, dirname, join as joinpath … … 24 20 HAVE_OPENCL = False 25 21 22 __all__ = [ 23 "list_models", "load_model_definition", "precompile_dll", 24 "load_model", "make_kernel", "call_kernel", "call_ER", "call_VR", 25 ] 26 26 27 27 def list_models(): -
sasmodels/data.py
rd15a908 r69ec80f 90 90 self.x, self.y, self.dx, self.dy = x, y, dx, dy 91 91 self.dxl = None 92 self.filename = None 93 self.qmin = x.min() if x is not None else np.NaN 94 self.qmax = x.max() if x is not None else np.NaN 95 self.mask = np.isnan(y) if y is not None else None 96 self._xaxis, self._xunit = "x", "" 97 self._yaxis, self._yunit = "y", "" 92 98 93 99 def xaxis(self, label, unit): … … 108 114 109 115 class Data2D(object): 110 def __init__(self): 116 def __init__(self, x=None, y=None, z=None, dx=None, dy=None, dz=None): 117 self.qx_data, self.dqx_data = x, dx 118 self.qy_data, self.dqy_data = y, dy 119 self.data, self.err_data = z, dz 120 self.mask = ~np.isnan(z) if z is not None else None 121 self.q_data = np.sqrt(x**2 + y**2) 122 self.qmin = 1e-16 123 self.qmax = np.inf 111 124 self.detector = [] 112 125 self.source = Source() 126 self.Q_unit = "1/A" 127 self.I_unit = "1/cm" 128 self.xaxis("Q_x", "A^{-1}") 129 self.yaxis("Q_y", "A^{-1}") 130 self.zaxis("Intensity", r"\text{cm}^{-1}") 131 self._xaxis, self._xunit = "x", "" 132 self._yaxis, self._yunit = "y", "" 133 self._zaxis, self._zunit = "z", "" 134 self.x_bins, self.y_bins = None, None 113 135 114 136 def xaxis(self, label, unit): … … 139 161 140 162 class Detector(object): 163 """ 164 Detector attributes. 165 """ 166 def __init__(self, pixel_size=(None, None), distance=None): 167 self.pixel_size = Vector(*pixel_size) 168 self.distance = distance 169 170 class Source(object): 171 """ 172 Beam attributes. 173 """ 141 174 def __init__(self): 142 self.pixel_size = Vector() 143 144 class Source(object): 145 pass 175 self.wavelength = np.NaN 176 self.wavelength_unit = "A" 146 177 147 178 … … 158 189 data = Data1D(q, Iq, dx=resolution * q, dy=dIq) 159 190 data.filename = "fake data" 160 data.qmin, data.qmax = q.min(), q.max()161 data.mask = np.zeros(len(q), dtype='bool')162 191 return data 163 192 … … 173 202 if qy is None: 174 203 qy = qx 204 # 5% dQ/Q resolution 175 205 Qx, Qy = np.meshgrid(qx, qy) 176 206 Qx, Qy = Qx.flatten(), Qy.flatten() 177 207 Iq = 100 * np.ones_like(Qx) 178 208 dIq = np.sqrt(Iq) 179 mask = np.ones(len(Iq), dtype='bool')180 181 data = Data2D()182 data.filename = "fake data"183 data.qx_data = Qx184 data.qy_data = Qy185 data.data = Iq186 data.err_data = dIq187 data.mask = mask188 data.qmin = 1e-16189 data.qmax = np.inf190 191 # 5% dQ/Q resolution192 209 if resolution != 0: 193 210 # https://www.ncnr.nist.gov/staff/hammouda/distance_learning/chapter_15.pdf … … 197 214 # radial (which instead it should be inverse). 198 215 Q = np.sqrt(Qx**2 + Qy**2) 199 d ata.dqx_data= resolution * Q200 d ata.dqy_data= resolution * Q216 dqx = resolution * Q 217 dqy = resolution * Q 201 218 else: 202 data.dqx_data = data.dqy_data = None 203 204 detector = Detector() 205 detector.pixel_size.x = 5 # mm 206 detector.pixel_size.y = 5 # mm 207 detector.distance = 4 # m 208 data.detector.append(detector) 219 dqx = dqy = None 220 221 data = Data2D(x=Qx, y=Qy, z=Iq, dx=dqx, dy=dqy, dz=dIq) 209 222 data.x_bins = qx 210 223 data.y_bins = qy 224 data.filename = "fake data" 225 226 # pixel_size in mm, distance in m 227 detector = Detector(pixel_size=(5, 5), distance=4) 228 data.detector.append(detector) 211 229 data.source.wavelength = 5 # angstroms 212 230 data.source.wavelength_unit = "A" 213 data.Q_unit = "1/A"214 data.I_unit = "1/cm"215 data.q_data = np.sqrt(Qx ** 2 + Qy ** 2)216 data.xaxis("Q_x", "A^{-1}")217 data.yaxis("Q_y", "A^{-1}")218 data.zaxis("Intensity", r"\text{cm}^{-1}")219 231 return data 220 232 … … 228 240 # do not repeat. 229 241 if hasattr(data, 'lam'): 230 _plot_result_sesans(data, None, None, plot_data=True, limits=limits)242 _plot_result_sesans(data, None, None, use_data=True, limits=limits) 231 243 elif hasattr(data, 'qx_data'): 232 _plot_result2D(data, None, None, view, plot_data=True, limits=limits)244 _plot_result2D(data, None, None, view, use_data=True, limits=limits) 233 245 else: 234 _plot_result1D(data, None, None, view, plot_data=True, limits=limits)246 _plot_result1D(data, None, None, view, use_data=True, limits=limits) 235 247 236 248 237 249 def plot_theory(data, theory, resid=None, view='log', 238 plot_data=True, limits=None):250 use_data=True, limits=None): 239 251 if hasattr(data, 'lam'): 240 _plot_result_sesans(data, theory, resid, plot_data=True, limits=limits)252 _plot_result_sesans(data, theory, resid, use_data=True, limits=limits) 241 253 elif hasattr(data, 'qx_data'): 242 _plot_result2D(data, theory, resid, view, plot_data, limits=limits)254 _plot_result2D(data, theory, resid, view, use_data, limits=limits) 243 255 else: 244 _plot_result1D(data, theory, resid, view, plot_data, limits=limits)256 _plot_result1D(data, theory, resid, view, use_data, limits=limits) 245 257 246 258 … … 251 263 except: 252 264 traceback.print_exc() 253 pass254 265 255 266 return wrapper … … 257 268 258 269 @protect 259 def _plot_result1D(data, theory, resid, view, plot_data, limits=None):270 def _plot_result1D(data, theory, resid, view, use_data, limits=None): 260 271 """ 261 272 Plot the data and residuals for 1D data. … … 264 275 from numpy.ma import masked_array, masked 265 276 266 plot_theory = theory is not None267 plot_resid = residis not None268 269 if data.y is None:270 plot_data = False 277 use_data = use_data and data.y is not None 278 use_theory = theory is not None 279 use_resid = resid is not None 280 num_plots = (use_data or use_theory) + use_resid 281 271 282 scale = data.x**4 if view == 'q4' else 1.0 272 283 273 if plot_data or plot_theory: 274 if plot_resid: 275 plt.subplot(121) 276 284 if use_data or use_theory: 277 285 #print(vmin, vmax) 278 286 all_positive = True 279 287 some_present = False 280 if plot_data:288 if use_data: 281 289 mdata = masked_array(data.y, data.mask.copy()) 282 290 mdata[~np.isfinite(mdata)] = masked … … 288 296 289 297 290 if plot_theory:298 if use_theory: 291 299 mtheory = masked_array(theory, data.mask.copy()) 292 300 mtheory[~np.isfinite(mtheory)] = masked … … 299 307 if limits is not None: 300 308 plt.ylim(*limits) 309 310 if num_plots > 1: 311 plt.subplot(1, num_plots, 1) 301 312 plt.xscale('linear' if not some_present else view) 302 313 plt.yscale('linear' … … 306 317 plt.ylabel('$I(q)$') 307 318 308 if plot_resid: 309 if plot_data or plot_theory: 310 plt.subplot(122) 311 319 if use_resid: 312 320 mresid = masked_array(resid, data.mask.copy()) 313 321 mresid[~np.isfinite(mresid)] = masked 314 322 some_present = (mresid.count() > 0) 323 324 if num_plots > 1: 325 plt.subplot(1, num_plots, (use_data or use_theory) + 1) 315 326 plt.plot(data.x/10, mresid, '-') 316 327 plt.xlabel("$q$/nm$^{-1}$") … … 320 331 321 332 @protect 322 def _plot_result_sesans(data, theory, resid, plot_data, limits=None):333 def _plot_result_sesans(data, theory, resid, use_data, limits=None): 323 334 import matplotlib.pyplot as plt 324 if data.y is None:325 plot_data = False326 plot_theory = theoryis not None327 plot_resid = resid is not None328 329 if plot_data or plot_theory:330 if plot_resid:331 plt.subplot(1 21)332 if plot_data:335 use_data = use_data and data.y is not None 336 use_theory = theory is not None 337 use_resid = resid is not None 338 num_plots = (use_data or use_theory) + use_resid 339 340 if use_data or use_theory: 341 if num_plots > 1: 342 plt.subplot(1, num_plots, 1) 343 if use_data: 333 344 plt.errorbar(data.x, data.y, yerr=data.dy) 334 345 if theory is not None: … … 340 351 341 352 if resid is not None: 342 if plot_data or plot_theory: 343 plt.subplot(122) 344 353 if num_plots > 1: 354 plt.subplot(1, num_plots, (use_data or use_theory) + 1) 345 355 plt.plot(data.x, resid, 'x') 346 356 plt.xlabel('spin echo length (nm)') … … 349 359 350 360 @protect 351 def _plot_result2D(data, theory, resid, view, plot_data, limits=None):361 def _plot_result2D(data, theory, resid, view, use_data, limits=None): 352 362 """ 353 363 Plot the data and residuals for 2D data. 354 364 """ 355 365 import matplotlib.pyplot as plt 356 if data.data is None:357 plot_data = False358 plot_theory = theoryis not None359 plot_resid = resid is not None366 use_data = use_data and data.data is not None 367 use_theory = theory is not None 368 use_resid = resid is not None 369 num_plots = use_data + use_theory + use_resid 360 370 361 371 # Put theory and data on a common colormap scale 362 if limits is None: 363 vmin, vmax = np.inf, -np.inf 364 if plot_data: 365 target = data.data[~data.mask] 366 datamin = target[target>0].min() if view == 'log' else target.min() 367 datamax = target.max() 368 vmin = min(vmin, datamin) 369 vmax = max(vmax, datamax) 370 if plot_theory: 371 theorymin = theory[theory>0].min() if view=='log' else theory.min() 372 theorymax = theory.max() 373 vmin = min(vmin, theorymin) 374 vmax = max(vmax, theorymax) 375 else: 372 vmin, vmax = np.inf, -np.inf 373 if use_data: 374 target = data.data[~data.mask] 375 datamin = target[target > 0].min() if view == 'log' else target.min() 376 datamax = target.max() 377 vmin = min(vmin, datamin) 378 vmax = max(vmax, datamax) 379 if use_theory: 380 theorymin = theory[theory > 0].min() if view == 'log' else theory.min() 381 theorymax = theory.max() 382 vmin = min(vmin, theorymin) 383 vmax = max(vmax, theorymax) 384 385 # Override data limits from the caller 386 if limits is not None: 376 387 vmin, vmax = limits 377 388 378 if plot_data: 379 if plot_theory and plot_resid: 380 plt.subplot(131) 381 elif plot_theory or plot_resid: 382 plt.subplot(121) 389 # Plot data 390 if use_data: 391 if num_plots > 1: 392 plt.subplot(1, num_plots, 1) 383 393 _plot_2d_signal(data, target, view=view, vmin=vmin, vmax=vmax) 384 394 plt.title('data') … … 386 396 h.set_label('$I(q)$') 387 397 388 if plot_theory: 389 if plot_data and plot_resid: 390 plt.subplot(132) 391 elif plot_data: 392 plt.subplot(122) 393 elif plot_resid: 394 plt.subplot(121) 398 # plot theory 399 if use_theory: 400 if num_plots > 1: 401 plt.subplot(1, num_plots, use_data+1) 395 402 _plot_2d_signal(data, theory, view=view, vmin=vmin, vmax=vmax) 396 403 plt.title('theory') … … 400 407 else '$I(q)$') 401 408 402 #if plot_data or plot_theory: 403 # plt.colorbar() 404 405 if plot_resid: 406 if plot_data and plot_theory: 407 plt.subplot(133) 408 elif plot_data or plot_theory: 409 plt.subplot(122) 409 # plot resid 410 if use_resid: 411 if num_plots > 1: 412 plt.subplot(1, num_plots, use_data+use_theory+1) 410 413 _plot_2d_signal(data, resid, view='linear') 411 414 plt.title('residuals') -
sasmodels/generate.py
re66c9f9 r190fc2b 197 197 # TODO: identify model files which have changed since loading and reload them. 198 198 199 __all__ = ["make", "doc", "sources", "convert_type"]200 201 199 import sys 202 200 from os.path import abspath, dirname, join as joinpath, exists, basename, \ … … 206 204 207 205 import numpy as np 206 207 __all__ = ["make", "doc", "sources", "convert_type"] 208 208 209 C_KERNEL_TEMPLATE_PATH = joinpath(dirname(__file__), 'kernel_template.c') 209 210 … … 216 217 F128 = None 217 218 218 219 219 # Scale and background, which are parameters common to every form factor 220 220 COMMON_PARAMETERS = [ … … 222 222 ["background", "1/cm", 0, [0, np.inf], "", "Source background"], 223 223 ] 224 225 224 226 225 # Conversion from units defined in the parameter table for each model … … 266 265 267 266 def format_units(units): 267 """ 268 Convert units into ReStructured Text format. 269 """ 268 270 return "string" if isinstance(units, list) else RST_UNITS.get(units, units) 269 271 … … 335 337 """ 336 338 Convert code from double precision to the desired type. 339 340 Floating point constants are tagged with 'f' for single precision or 'L' 341 for long double precision. 337 342 """ 338 343 if dtype == F16: … … 350 355 351 356 def _convert_type(source, type_name, constant_flag): 357 """ 358 Replace 'double' with *type_name* in *source*, tagging floating point 359 constants with *constant_flag*. 360 """ 352 361 # Convert double keyword to float/long double/half. 353 362 # Accept an 'n' # parameter for vector # values, where n is 2, 4, 8 or 16. … … 625 634 %re.escape(string.punctuation)) 626 635 def _convert_section_titles_to_boldface(lines): 636 """ 637 Do the actual work of identifying and converting section headings. 638 """ 627 639 prior = None 628 640 for line in lines: … … 642 654 yield prior 643 655 644 def convert_section_titles_to_boldface(string): 645 return "\n".join(_convert_section_titles_to_boldface(string.split('\n'))) 656 def convert_section_titles_to_boldface(s): 657 """ 658 Use explicit bold-face rather than section headings so that the table of 659 contents is not polluted with section names from the model documentation. 660 661 Sections are identified as the title line followed by a line of punctuation 662 at least as long as the title line. 663 """ 664 return "\n".join(_convert_section_titles_to_boldface(s.split('\n'))) 646 665 647 666 def doc(kernel_module): … … 666 685 667 686 def demo_time(): 687 """ 688 Show how long it takes to process a model. 689 """ 668 690 from .models import cylinder 669 691 import datetime … … 674 696 675 697 def main(): 698 """ 699 Program which prints the source produced by the model. 700 """ 676 701 if len(sys.argv) <= 1: 677 702 print("usage: python -m sasmodels.generate modelname") -
sasmodels/models/lamellarCailleHG.py
rd18f8a8 rd4666ca 92 92 93 93 parameters = [ 94 95 96 [ "tail_length", "Ang",10, [0, inf], "volume",97 "Tail thickness"],98 [ "head_length", "Ang",2, [0, inf], "volume",99 "head thickness"],100 [ "Nlayers", "",30, [0, inf], "",101 "Number of layers"],102 [ "spacing", "Ang", 40., [0.0,inf], "volume",103 "d-spacing of Caille S(Q)"],104 [ "Caille_parameter", "", 0.001, [0.0,0.8], "",105 "Caille parameter"],106 [ "sld", "1e-6/Ang^2", 0.4, [-inf,inf], "",107 "Tail scattering length density"],108 [ "head_sld", "1e-6/Ang^2", 2.0, [-inf,inf], "",109 "Head scattering length density"],110 [ "solvent_sld", "1e-6/Ang^2", 6, [-inf,inf], "",111 "Solvent scattering length density"],94 # [ "name", "units", default, [lower, upper], "type", 95 # "description" ], 96 ["tail_length", "Ang", 10, [0, inf], "volume", 97 "Tail thickness"], 98 ["head_length", "Ang", 2, [0, inf], "volume", 99 "head thickness"], 100 ["Nlayers", "", 30, [0, inf], "", 101 "Number of layers"], 102 ["spacing", "Ang", 40., [0.0, inf], "volume", 103 "d-spacing of Caille S(Q)"], 104 ["Caille_parameter", "", 0.001, [0.0, 0.8], "", 105 "Caille parameter"], 106 ["sld", "1e-6/Ang^2", 0.4, [-inf, inf], "", 107 "Tail scattering length density"], 108 ["head_sld", "1e-6/Ang^2", 2.0, [-inf, inf], "", 109 "Head scattering length density"], 110 ["solvent_sld", "1e-6/Ang^2", 6, [-inf, inf], "", 111 "Solvent scattering length density"], 112 112 ] 113 113 114 source = [ 114 source = ["lamellarCailleHG_kernel.c"] 115 115 116 116 # No volume normalization despite having a volume parameter … … 128 128 129 129 demo = dict( 130 scale=1, background=0, 131 Nlayers=20, 132 spacing=200., Caille_parameter=0.05, 133 tail_length=15,head_length=10, 134 #sld=-1, head_sld=4.0, solvent_sld=6.0, 135 sld=-1, head_sld=4.1, solvent_sld=6.0, 136 tail_length_pd= 0.1, tail_length_pd_n=20, 137 head_length_pd= 0.05, head_length_pd_n=30, 138 spacing_pd= 0.2, spacing_pd_n=40 139 ) 130 scale=1, background=0, 131 Nlayers=20, spacing=200., Caille_parameter=0.05, 132 tail_length=15, head_length=10, 133 #sld=-1, head_sld=4.0, solvent_sld=6.0, 134 sld=-1, head_sld=4.1, solvent_sld=6.0, 135 tail_length_pd=0.1, tail_length_pd_n=20, 136 head_length_pd=0.05, head_length_pd_n=30, 137 spacing_pd=0.2, spacing_pd_n=40, 138 ) 140 139 141 140 oldname = 'LamellarPSHGModel' 142 oldpars = dict(tail_length='deltaT',head_length='deltaH',Nlayers='n_plates',Caille_parameter='caille', sld='sld_tail', head_sld='sld_head',solvent_sld='sld_solvent') 141 oldpars = dict( 142 tail_length='deltaT', head_length='deltaH', Nlayers='n_plates', 143 Caille_parameter='caille', sld='sld_tail', head_sld='sld_head', 144 solvent_sld='sld_solvent') -
sasmodels/models/lib/sph_j1c.c
r9c461c7 ra3e78c3 1 1 /** 2 * Spherical Bessel function j1(x)/x2 * Spherical Bessel function 3*j1(x)/x 3 3 * 4 4 * Used for low q to avoid cancellation error. … … 8 8 * requires the first 3 terms. Double precision requires the 4th term. 9 9 * The fifth term is not needed, and is commented out. 10 * Taylor expansion: 11 * 1.0 + q2*(-3./30. + q2*(3./840.))+ q2*(-3./45360. + q2*(3./3991680.)))) 12 * Expression returned from Herbie (herbie.uwpise.org/demo): 13 * const double t = ((1. + 3.*q2*q2/5600.) - q2/20.); 14 * return t*t; 10 15 */ 11 16 … … 18 23 SINCOS(q, sin_q, cos_q); 19 24 20 const double bessel = (q < 1.e-1) ?21 1.0 + q2*(-3./30. + q2*(3./840. + q2*(-3./45360.))) // + q2*(3./3991680.)))25 const double bessel = (q < 0.384038453352533) 26 ? (1.0 + q2*(-3./30. + q2*(3./840.))) 22 27 : 3.0*(sin_q/q - cos_q)/q2; 23 28 24 29 return bessel; 30 31 /* 32 // Code to test various expressions 33 if (sizeof(q2) > 4) { 34 return 3.0*(sin_q/q - cos_q)/q2; 35 } else if (q < 0.384038453352533) { 36 //const double t = ((1. + 3.*q2*q2/5600.) - q2/20.); return t*t; 37 return 1.0 + q2*q2*(3./840.) - q2*(3./30.); 38 //return 1.0 + q2*(-3./30. + q2*(3./840.)); 39 //return 1.0 + q2*(-3./30. + q2*(3./840. + q2*(-3./45360.))); 40 //return 1.0 + q2*(-3./30. + q2*(3./840. + q2*(-3./45360. + q2*(3./3991680.)))); 41 } else { 42 return 3.0*(sin_q/q - cos_q)/q2; 43 } 44 */ 25 45 } -
sasmodels/models/pearl_necklace.py
rcf404cb r841753c 1 1 r""" 2 This model provides the form factor for a pearl necklace composed of two 3 elements: *N* pearls (homogeneous spheres of radius *R*) freely jointed by *M* 2 This model provides the form factor for a pearl necklace composed of two 3 elements: *N* pearls (homogeneous spheres of radius *R*) freely jointed by *M* 4 4 rods (like strings - with a total mass *Mw* = *M* \* *m*\ :sub:`r` + *N* \* *m*\ 5 :sub:`s`, and the string segment length (or edge separation) *l* 5 :sub:`s`, and the string segment length (or edge separation) *l* 6 6 (= *A* - 2\ *R*)). *A* is the center-to-center pearl separation distance. 7 7 … … 13 13 ---------- 14 14 15 The output of the scattering intensity function for the PearlNecklaceModel is 15 The output of the scattering intensity function for the PearlNecklaceModel is 16 16 given by (Schweins, 2004) 17 17 … … 37 37 \beta(q) &= \frac{\int_{qR}^{q(A-R)}\frac{sin(t)}{t}dt}{ql} 38 38 39 where the mass *m*\ :sub:`i` is (SLD\ :sub:`i` - SLD\ :sub:`solvent`) \* 39 where the mass *m*\ :sub:`i` is (SLD\ :sub:`i` - SLD\ :sub:`solvent`) \* 40 40 (volume of the *N* pearls/rods). *V* is the total volume of the necklace. 41 41 42 The 2D scattering intensity is the same as $P(q)$ above, regardless of the 42 The 2D scattering intensity is the same as $P(q)$ above, regardless of the 43 43 orientation of the *q* vector. 44 44 … … 54 54 REFERENCE 55 55 56 R Schweins and K Huber, *Particle Scattering Factor of Pearl Necklace Chains*, 56 R Schweins and K Huber, *Particle Scattering Factor of Pearl Necklace Chains*, 57 57 *Macromol. Symp.* 211 (2004) 25-42 2004 58 58 """ … … 79 79 80 80 # ["name", "units", default, [lower, upper], "type","description"], 81 parameters = [["radius", "Angstrom", 80.0, [0, inf], "volume", 81 parameters = [["radius", "Angstrom", 80.0, [0, inf], "volume", 82 82 "Mean radius of the chained spheres"], 83 ["edge_separation", "Angstrom", 350.0, [0, inf], "volume", 83 ["edge_separation", "Angstrom", 350.0, [0, inf], "volume", 84 84 "Mean separation of chained particles"], 85 ["string_thickness", "Angstrom", 2.5, [0, inf], "volume", 85 ["string_thickness", "Angstrom", 2.5, [0, inf], "volume", 86 86 "Thickness of the chain linkage"], 87 ["number_of_pearls", "none", 3, [0, inf], "volume", 87 ["number_of_pearls", "none", 3, [0, inf], "volume", 88 88 "Mean number of pearls in each necklace"], 89 ["sld", "Angstrom^2", 1.0, [-inf, inf], "", 89 ["sld", "Angstrom^2", 1.0, [-inf, inf], "", 90 90 "Scattering length density of the chained spheres"], 91 ["string_sld", "Angstrom^2", 1.0, [-inf, inf], "", 91 ["string_sld", "Angstrom^2", 1.0, [-inf, inf], "", 92 92 "Scattering length density of the chain linkage"], 93 ["solvent_sld", "Angstrom^2", 6.3, [-inf, inf], "", 93 ["solvent_sld", "Angstrom^2", 6.3, [-inf, inf], "", 94 94 "Scattering length density of the solvent"], 95 95 ] 96 96 97 97 source = ["lib/Si.c", "pearl_necklace.c"] … … 110 110 # names and the target sasview model name. 111 111 oldname = 'PearlNecklaceModel' 112 oldpars = dict(scale='scale', background='background',radius='radius',112 oldpars = dict(scale='scale', background='background', radius='radius', 113 113 number_of_pearls='num_pearls', solvent_sld='sld_solv', 114 114 string_thickness='thick_string', sld='sld_pearl', -
sasmodels/models/power_law.py
reb69cce r841753c 12 12 I(q) = \text{scale} \cdot q^{-\text{power}} + \text{background} 13 13 14 Note the minus sign in front of the exponent. The exponent *power* 14 Note the minus sign in front of the exponent. The exponent *power* 15 15 should therefore be entered as a **positive** number for fitting. 16 16 17 Also note that unlike many other models, *scale* in this model 18 is NOT explicitly related to a volume fraction. Be careful if 17 Also note that unlike many other models, *scale* in this model 18 is NOT explicitly related to a volume fraction. Be careful if 19 19 combining this model with other models. 20 20 … … 31 31 from numpy import inf, sqrt 32 32 33 name = 33 name = "power_law" 34 34 title = "Simple power law with a flat background" 35 35 36 description = """ \37 38 39 40 36 description = """ 37 Evaluates the function 38 I(q) = scale * q^(-power) + background 39 NB: enter power as a positive number! 40 """ 41 41 category = "shape-independent" 42 42 … … 45 45 46 46 # NB: Scale and Background are implicit parameters on every model 47 def Iq(q,power): 47 def Iq(q, power): 48 # pylint: disable=missing-docstring 48 49 inten = (q**-power) 49 50 return inten … … 51 52 52 53 def Iqxy(qx, qy, *args): 54 # pylint: disable=missing-docstring 53 55 return Iq(sqrt(qx ** 2 + qy ** 2), *args) 54 56 Iqxy.vectorized = True # Iqxy accepts an array of qx, qy values … … 64 66 65 67 tests = [ 66 [ {'scale': 1.0, 'power': 4.0, 'background' : 0.0}, [0.0106939, 0.469418], [7.64644e+07, 20.5949]] 67 ] 68 [{'scale': 1.0, 'power': 4.0, 'background' : 0.0}, 69 [0.0106939, 0.469418], [7.64644e+07, 20.5949]], 70 ] -
sasmodels/resolution.py
rfdc538a r190fc2b 5 5 """ 6 6 from __future__ import division 7 8 from scipy.special import erf 9 from numpy import sqrt, log, log10 10 import numpy as np 7 11 8 12 __all__ = ["Resolution", "Perfect1D", "Pinhole1D", "Slit1D", … … 11 15 "interpolate", "linear_extrapolation", "geometric_extrapolation", 12 16 ] 13 14 from scipy.special import erf15 from numpy import sqrt, log, log1016 import numpy as np17 17 18 18 MINIMUM_RESOLUTION = 1e-8 -
sasmodels/resolution2d.py
r9404dd3 r841753c 2 2 #This software was developed by the University of Tennessee as part of the 3 3 #Distributed Data Analysis of Neutron Scattering Experiments (DANSE) 4 #project funded by the US National Science Foundation. 4 #project funded by the US National Science Foundation. 5 5 #See the license text in license.txt 6 6 """ … … 19 19 ## Defaults 20 20 NR = {'xhigh':10, 'high':5, 'med':5, 'low':3} 21 NPHI = {'xhigh':20, 'high':12, 'med':6, 'low':4}21 NPHI = {'xhigh':20, 'high':12, 'med':6, 'low':4} 22 22 23 23 class Pinhole2D(Resolution): … … 25 25 Gaussian Q smearing class for SAS 2d data 26 26 """ 27 27 28 28 def __init__(self, data=None, index=None, 29 29 nsigma=NSIGMA, accuracy='Low', coords='polar'): 30 30 """ 31 31 Assumption: equally spaced bins in dq_r, dq_phi space. 32 32 33 33 :param data: 2d data used to set the smearing parameters 34 34 :param index: 1d array with len(data) to define the range … … 42 42 ## number of bins in r axis for over-sampling 43 43 self.nr = NR[accuracy.lower()] 44 ## number of bins in phi axis for over-sampling 44 ## number of bins in phi axis for over-sampling 45 45 self.nphi = NPHI[accuracy.lower()] 46 46 ## maximum nsigmas 47 self.nsigma = nsigma47 self.nsigma = nsigma 48 48 self.coords = coords 49 49 self._init_data(data, index) … … 85 85 def _calc_res(self): 86 86 """ 87 Over sampling of r_nbins times phi_nbins, calculate Gaussian weights, 87 Over sampling of r_nbins times phi_nbins, calculate Gaussian weights, 88 88 then find smeared intensity 89 """ 89 """ 90 90 nr, nphi = self.nr, self.nphi 91 91 # Total number of bins = # of bins … … 143 143 if self.coords == 'polar': 144 144 q_r = sqrt(qx**2 + qy**2) 145 qx_res = ( 145 qx_res = ((dqx*cos(dphi) + q_r) * cos(-q_phi) + 146 146 dqy*sin(dphi) * sin(-q_phi)) 147 147 qy_res = (-(dqx*cos(dphi) + q_r) * sin(-q_phi) + … … 167 167 else: 168 168 return theory 169 170 """171 if __name__ == '__main__':172 ## Test w/ 2D linear function173 x = 0.001*np.arange(1, 11)174 dx = np.ones(len(x))*0.0003175 y = 0.001*np.arange(1, 11)176 dy = np.ones(len(x))*0.001177 z = np.ones(10)178 dz = sqrt(z)179 180 from sas.dataloader import Data2D181 #for i in range(10): print(i, 0.001 + i*0.008/9.0)182 #for i in range(100): print(i, int(math.floor( (i/ (100/9.0)) )))183 out = Data2D()184 out.data = z185 out.qx_data = x186 out.qy_data = y187 out.dqx_data = dx188 out.dqy_data = dy189 out.q_data = sqrt(dx * dx + dy * dy)190 index = np.ones(len(x), dtype = bool)191 out.mask = index192 from sas.models.LineModel import LineModel193 model = LineModel()194 model.setParam("A", 0)195 196 smear = Smearer2D(out, model, index)197 #smear.set_accuracy('Xhigh')198 value = smear.get_value()199 ## All data are ones, so the smeared should also be ones.200 print("Data length =", len(value))201 print(" 2D linear function, I = 0 + 1*qy")202 text = " Gaussian weighted averaging on a 2D linear function will "203 text += "provides the results same as without the averaging."204 print(text)205 print("qx_data", "qy_data", "I_nonsmear", "I_smeared")206 for ind in range(len(value)):207 print(x[ind], y[ind], model.evalDistribution([x, y])[ind], value[ind])208 209 210 if __name__ == '__main__':211 ## Another Test w/ constant function212 x = 0.001*np.arange(1,11)213 dx = np.ones(len(x))*0.001214 y = 0.001*np.arange(1,11)215 dy = np.ones(len(x))*0.001216 z = np.ones(10)217 dz = sqrt(z)218 219 from DataLoader import Data2D220 #for i in range(10): print(i, 0.001 + i*0.008/9.0)221 #for i in range(100): print(i, int(math.floor( (i/ (100/9.0)) )))222 out = Data2D()223 out.data = z224 out.qx_data = x225 out.qy_data = y226 out.dqx_data = dx227 out.dqy_data = dy228 index = np.ones(len(x), dtype = bool)229 out.mask = index230 from sas.models.Constant import Constant231 model = Constant()232 233 value = Smearer2D(out,model,index).get_value()234 ## All data are ones, so the smeared values should also be ones.235 print("Data length =",len(value), ", Data=",value)236 """ -
sasmodels/sesans.py
r384d114 r190fc2b 43 43 *I* [cm$^{-1}$] is the value of the SANS model at *q* 44 44 """ 45 G = np.zeros (len(SElength), 'd')46 for i in range(len(SElength)):47 integr = besselj(0, q*SElength[i])*Iq*q48 G[i] = np.sum(integr )45 G = np.zeros_like(SElength, 'd') 46 for i, SElength_i in enumerate(SElength): 47 integral = besselj(0, q*SElength_i)*Iq*q 48 G[i] = np.sum(integral) 49 49 50 50 # [m^-1] step size in q, needed for integration 51 dq =(q[1]-q[0])*1e1051 dq = (q[1]-q[0])*1e10 52 52 53 53 # integration step, convert q into [m**-1] and 2 pi circle integration -
sasmodels/models/be_polyelectrolyte.py
r841753c r168052c 44 44 J F Joanny, L Leibler, *Journal de Physique*, 51 (1990) 545 45 45 46 A Moussaid, F Schosseler, J P Munch, S Candau, *J. Journal de Physique II France*, 3 (1993) 573 46 A Moussaid, F Schosseler, J P Munch, S Candau, 47 *J. Journal de Physique II France*, 3 (1993) 573 47 48 48 49 E Raphael, J F Joanny, *Europhysics Letters*, 11 (1990) 179 … … 55 56 title = "Polyelectrolyte with the RPA expression derived by Borue and Erukhimovich" 56 57 description = """ 57 Evaluate58 F(x) = K 1/(4 pi Lb (alpha)^(2)) (q^(2)+k2)/(1+(r02)^(2))59 (q^(2)+k2) (q^(2)-(12 h C/b^(2)))58 Evaluate 59 F(x) = K 1/(4 pi Lb (alpha)^(2)) (q^(2)+k2)/(1+(r02)^(2)) 60 (q^(2)+k2) (q^(2)-(12 h C/b^(2))) 60 61 61 has 3 internal parameters :62 The inverse Debye Length: K2 = 4 pi Lb (2 Cs+alpha C)63 r02 =1/alpha/Ca^(0.5) (B/(48 pi Lb)^(0.5))64 Ca = 6.022136e-4 C65 """62 has 3 internal parameters : 63 The inverse Debye Length: K2 = 4 pi Lb (2 Cs+alpha C) 64 r02 =1/alpha/Ca^(0.5) (B/(48 pi Lb)^(0.5)) 65 Ca = 6.022136e-4 C 66 """ 66 67 category = "shape-independent" 67 68 68 # pylint: disable=bad-whitespace, line-too-long69 # ["name", 69 # pylint: disable=bad-whitespace, line-too-long 70 # ["name", "units", default, [lower, upper], "type", "description"], 70 71 parameters = [ 71 ["contrast_factor", "barns", 72 ["bjerrum_length", "Ang", 73 ["virial_param", "1/Ang^2", 74 ["monomer_length", "Ang", 75 ["salt_concentration", "mol/L", 76 ["ionization_degree", "", 77 ["polymer_concentration", "mol/L", 72 ["contrast_factor", "barns", 10.0, [-inf, inf], "", "Contrast factor of the polymer"], 73 ["bjerrum_length", "Ang", 7.1, [0, inf], "", "Bjerrum length"], 74 ["virial_param", "1/Ang^2", 12.0, [-inf, inf], "", "Virial parameter"], 75 ["monomer_length", "Ang", 10.0, [0, inf], "", "Monomer length"], 76 ["salt_concentration", "mol/L", 0.0, [-inf, inf], "", "Concentration of monovalent salt"], 77 ["ionization_degree", "", 0.05, [0, inf], "", "Degree of ionization"], 78 ["polymer_concentration", "mol/L", 0.7, [0, inf], "", "Polymer molar concentration"], 78 79 ] 79 # pylint: enable=bad-whitespace, line-too-long80 # pylint: enable=bad-whitespace, line-too-long 80 81 81 82 82 83 def Iq(q, 83 contrast_factor, 84 bjerrum_length, 85 virial_param, 86 monomer_length, 87 salt_concentration, 88 ionization_degree, 89 polymer_concentration): 84 contrast_factor=10.0, 85 bjerrum_length=7.1, 86 virial_param=12.0, 87 monomer_length=10.0, 88 salt_concentration=0.0, 89 ionization_degree=0.05, 90 polymer_concentration=0.7): 91 """ 92 :param q: Input q-value 93 :param contrast_factor: Contrast factor of the polymer 94 :param bjerrum_length: Bjerrum length 95 :param virial_param: Virial parameter 96 :param monomer_length: Monomer length 97 :param salt_concentration: Concentration of monovalent salt 98 :param ionization_degree: Degree of ionization 99 :param polymer_concentration: Polymer molar concentration 100 :return: 1-D intensity 101 """ 90 102 91 103 concentration = polymer_concentration * 6.022136e-4 92 104 93 105 k_square = 4.0 * pi * bjerrum_length * (2*salt_concentration + 94 ionization_degree * concentration)106 ionization_degree * concentration) 95 107 96 108 r0_square = 1.0/ionization_degree/sqrt(concentration) * \ 97 (monomer_length/sqrt((48.0*pi*bjerrum_length)))109 (monomer_length/sqrt((48.0*pi*bjerrum_length))) 98 110 99 111 term1 = contrast_factor/(4.0 * pi * bjerrum_length * 100 ionization_degree**2) * (q**2 + k_square)112 ionization_degree**2) * (q**2 + k_square) 101 113 102 114 term2 = 1.0 + r0_square**2 * (q**2 + k_square) * \ … … 109 121 110 122 def Iqxy(qx, qy, *args): 111 return Iq(sqrt(qx**2 + qy**2), *args) 123 """ 124 :param qx: Input q_x-value 125 :param qy: Input q_y-value 126 :param args: Remaining arguments 127 :return: 2D-Intensity 128 """ 129 iq = Iq(sqrt(qx**2 + qy**2), *args) 130 return iq 112 131 113 132 Iqxy.vectorized = True # Iqxy accepts an array of qx, qy values … … 134 153 polymer_concentration='c') 135 154 136 # pylint: disable=bad-whitespace137 155 tests = [ 156 138 157 # Accuracy tests based on content in test/utest_other_models.py 139 158 [{'contrast_factor': 10.0, … … 176 195 }, 200., 1.80664667511e-06], 177 196 ] 178 # pylint: enable=bad-whitespace -
sasmodels/models/flexible_cylinder.py
rf94d8a2 r168052c 1 1 r""" 2 This model provides the form factor, $P(q)$, for a flexible cylinder where the form factor3 is normalized by the volume of the cylinder.2 This model provides the form factor, $P(q)$, for a flexible cylinder 3 where the form factor is normalized by the volume of the cylinder. 4 4 **Inter-cylinder interactions are NOT provided for.** 5 5 … … 8 8 P(q) = \text{scale} \left<F^2\right>/V + \text{background} 9 9 10 where the averaging $\left<\ldots\right>$ is applied only for the 1D calculation 10 where the averaging $\left<\ldots\right>$ is applied only for the 1D 11 calculation 11 12 12 The 2D scattering intensity is the same as 1D, regardless of the orientation of the q vector which is defined as 13 The 2D scattering intensity is the same as 1D, regardless of the orientation of 14 the q vector which is defined as 13 15 14 16 .. math:: … … 22 24 23 25 24 The chain of contour length, $L$, (the total length) can be described as a chain of some number of 25 locally stiff segments of length $l_p$, the persistence length (the length along the cylinder over 26 which the flexible cylinder can be considered a rigid rod). 26 The chain of contour length, $L$, (the total length) can be described as a 27 chain of some number of locally stiff segments of length $l_p$, the persistence 28 length (the length along the cylinder over which the flexible cylinder can be 29 considered a rigid rod). 27 30 The Kuhn length $(b = 2*l_p)$ is also used to describe the stiffness of a chain. 28 31 29 32 The returned value is in units of $cm^-1$, on absolute scale. 30 33 31 In the parameters, the sldCyl and sldSolv represent the SLD of the chain/cylinder and solvent respectively. 34 In the parameters, the sldCyl and sldSolv represent the SLD of the chain/cylinder 35 and solvent respectively. 32 36 33 37 … … 37 41 38 42 39 Our model uses the form factor calculations implemented in a c-library provided by the NIST Center40 for Neutron Research (Kline, 2006).43 Our model uses the form factor calculations implemented in a c-library provided 44 by the NIST Center for Neutron Research (Kline, 2006). 41 45 42 46 … … 45 49 'Method 3 With Excluded Volume' is used. 46 50 The model is a parametrization of simulations of a discrete representation 47 of the worm-like chain model of Kratky and Porod applied in the pseudocontinuous limit. 51 of the worm-like chain model of Kratky and Porod applied in the 52 pseudocontinuous limit. 48 53 See equations (13,26-27) in the original reference for the details. 49 54 … … 51 56 ---------- 52 57 53 J S Pedersen and P Schurtenberger. *Scattering functions of semiflexible polymers with and 54 without excluded volume effects.* Macromolecules, 29 (1996) 7602-7612 58 J S Pedersen and P Schurtenberger. *Scattering functions of semiflexible 59 polymers with and without excluded volume effects.* Macromolecules, 60 29 (1996) 7602-7612 55 61 56 62 Correction of the formula can be found in 57 63 58 W R Chen, P D Butler and L J Magid, *Incorporating Intermicellar Interactions in the Fitting of 59 SANS Data from Cationic Wormlike Micelles.* Langmuir, 22(15) 2006 6539-6548 64 W R Chen, P D Butler and L J Magid, *Incorporating Intermicellar Interactions 65 in the Fitting of SANS Data from Cationic Wormlike Micelles.* Langmuir, 66 22(15) 2006 6539-6548 60 67 """ 61 68 from numpy import inf 62 69 63 70 name = "flexible_cylinder" 64 title = "Flexible cylinder where the form factor is normalized by the volume of the cylinder." 65 description = """Note : scale and contrast=sld-solvent_sld are both multiplicative factors in the 66 model and are perfectly correlated. One or 67 both of these parameters must be held fixed 71 title = "Flexible cylinder where the form factor is normalized by the volume" \ 72 "of the cylinder." 73 description = """Note : scale and contrast=sld-solvent_sld are both 74 multiplicative factors in the model and are perfectly 75 correlated. One or both of these parameters must be held fixed 68 76 during model fitting. 69 77 """ … … 71 79 category = "shape:cylinder" 72 80 81 # pylint: disable=bad-whitespace, line-too-long 73 82 # ["name", "units", default, [lower, upper], "type", "description"], 74 83 parameters = [ 75 76 77 78 79 80 81 84 ["length", "Ang", 1000.0, [0, inf], "volume", "Length of the flexible cylinder"], 85 ["kuhn_length", "Ang", 100.0, [0, inf], "volume", "Kuhn length of the flexible cylinder"], 86 ["radius", "Ang", 20.0, [0, inf], "volume", "Radius of the flexible cylinder"], 87 ["sld", "1e-6/Ang^2", 1.0, [-inf, inf], "", "Cylinder scattering length density"], 88 ["solvent_sld", "1e-6/Ang^2", 6.3, [-inf, inf], "", "Solvent scattering length density"], 89 ] 90 # pylint: enable=bad-whitespace, line-too-long 82 91 source = ["lib/J1.c", "lib/wrc_cyl.c", "flexible_cylinder.c"] 83 92 … … 94 103 95 104 tests = [ 96 97 98 99 100 101 102 103 104 105 # Accuracy tests based on content in test/utest_other_models.py 106 # Currently fails in OCL 107 # [{'length': 1000.0, 108 # 'kuhn_length': 100.0, 109 # 'radius': 20.0, 110 # 'sld': 1.0, 111 # 'solvent_sld': 6.3, 112 # 'background': 0.0001, 113 # }, 0.001, 3509.2187], 105 114 106 107 [{'length':1000.0,108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 115 # Additional tests with larger range of parameters 116 [{'length': 1000.0, 117 'kuhn_length': 100.0, 118 'radius': 20.0, 119 'sld': 1.0, 120 'solvent_sld': 6.3, 121 'background': 0.0001, 122 }, 1.0, 0.000595345], 123 [{'length': 10.0, 124 'kuhn_length': 800.0, 125 'radius': 2.0, 126 'sld': 6.0, 127 'solvent_sld': 12.3, 128 'background': 0.001, 129 }, 0.1, 1.55228], 130 [{'length': 100.0, 131 'kuhn_length': 800.0, 132 'radius': 50.0, 133 'sld': 0.1, 134 'solvent_sld': 5.1, 135 'background': 0.0, 136 }, 1.0, 0.000938456] 137 ] 129 138 -
sasmodels/models/gauss_lorentz_gel.py
r07a6700 r168052c 42 42 ---------- 43 43 44 G Evmenenko, E Theunissen, K Mortensen, H Reynaers, *Polymer*, 42 (2001) 2907-2913 44 G Evmenenko, E Theunissen, K Mortensen, H Reynaers, *Polymer*, 45 42 (2001) 2907-2913 45 46 46 47 """ 47 48 48 from numpy import inf, pi,sqrt, exp49 from numpy import inf, sqrt, exp 49 50 50 51 name = "gauss_lorentz_gel" … … 63 64 """ 64 65 category = "shape-independent" 65 66 # pylint: disable=bad-whitespace, line-too-long 66 67 # ["name", "units", default, [lower, upper], "type", "description"], 67 68 parameters = [["gauss_scale_factor", "", 100.0, [-inf, inf], "", "Gauss scale factor"], … … 69 70 ["lorentz_scale_factor", "", 50.0, [-inf, inf], "", "Lorentzian scale factor"], 70 71 ["dynamic_cor_length", "Ang", 20.0, [0, inf], "", "Dynamic correlation length"], 71 72 72 ] 73 # pylint: enable=bad-whitespace, line-too-long 73 74 74 75 def Iq(q, 75 gauss_scale_factor, 76 static_cor_length, 77 lorentz_scale_factor, 78 dynamic_cor_length): 76 gauss_scale_factor=100.0, 77 static_cor_length=100.0, 78 lorentz_scale_factor=50.0, 79 dynamic_cor_length=20.0): 80 """ 79 81 80 term1 = gauss_scale_factor *\ 81 exp(-1.0*q*q*static_cor_length*static_cor_length/2.0) 82 term2 = lorentz_scale_factor /\ 83 (1.0+(q*dynamic_cor_length)*(q*dynamic_cor_length)) 82 :param q: Input q-value 83 :param gauss_scale_factor: Gauss scale factor 84 :param static_cor_length: Static correlation length 85 :param lorentz_scale_factor: Lorentzian scale factor 86 :param dynamic_cor_length: Dynamic correlation length 87 :return: 1-D intensity 88 """ 84 89 85 return term1 + term2 90 term1 = gauss_scale_factor *\ 91 exp(-1.0*q*q*static_cor_length*static_cor_length/2.0) 92 term2 = lorentz_scale_factor /\ 93 (1.0+(q*dynamic_cor_length)*(q*dynamic_cor_length)) 94 95 return term1 + term2 86 96 87 97 Iq.vectorized = True # Iq accepts an array of q values … … 89 99 90 100 def Iqxy(qx, qy, *args): 91 iq = Iq(sqrt(qx**2 + qy**2), *args) 101 """ 102 :param qx: Input q_x-value 103 :param qy: Input q_y-value 104 :param args: Remaining aruments 105 :return: 2-D intensity 106 """ 92 107 93 return iq108 return Iq(sqrt(qx**2 + qy**2), *args) 94 109 95 110 Iqxy.vectorized = True # Iqxy accepts an array of qx, qy values … … 111 126 112 127 tests = [ 113 # Accuracy tests based on content in test/utest_extra_models.py114 [{'gauss_scale_factor': 100.0,115 'static_cor_length': 100.0,116 'lorentz_scale_factor': 50.0,117 'dynamic_cor_length': 20.0,118 }, 0.001, 149.481],119 128 120 [{'gauss_scale_factor': 100.0, 121 'static_cor_length': 100.0, 122 'lorentz_scale_factor': 50.0, 123 'dynamic_cor_length': 20.0, 124 }, 0.105363, 9.1903], 129 # Accuracy tests based on content in test/utest_extra_models.py 130 [{'gauss_scale_factor': 100.0, 131 'static_cor_length': 100.0, 132 'lorentz_scale_factor': 50.0, 133 'dynamic_cor_length': 20.0, 134 }, 0.001, 149.481], 125 135 126 127 128 129 130 }, 0.441623, 0.632811],136 [{'gauss_scale_factor': 100.0, 137 'static_cor_length': 100.0, 138 'lorentz_scale_factor': 50.0, 139 'dynamic_cor_length': 20.0, 140 }, 0.105363, 9.1903], 131 141 132 # Additional tests with larger range of parameters 133 [{'gauss_scale_factor': 10.0, 134 'static_cor_length': 100.0, 135 'lorentz_scale_factor': 3.0, 136 'dynamic_cor_length': 1.0, 137 }, 0.1, 2.9702970297], 142 [{'gauss_scale_factor': 100.0, 143 'static_cor_length': 100.0, 144 'lorentz_scale_factor': 50.0, 145 'dynamic_cor_length': 20.0, 146 }, 0.441623, 0.632811], 138 147 139 [{'gauss_scale_factor': 10.0,140 'static_cor_length': 100.0,141 'lorentz_scale_factor': 3.0,142 'dynamic_cor_length': 1.0,143 'background': 100.0144 }, 5.0, 100.115384615],148 # Additional tests with larger range of parameters 149 [{'gauss_scale_factor': 10.0, 150 'static_cor_length': 100.0, 151 'lorentz_scale_factor': 3.0, 152 'dynamic_cor_length': 1.0, 153 }, 0.1, 2.9702970297], 145 154 146 [{'gauss_scale_factor': 10.0, 147 'static_cor_length': 100.0, 148 'lorentz_scale_factor': 3.0, 149 'dynamic_cor_length': 1.0, 150 }, 200., 7.49981250469e-05], 151 ] 155 [{'gauss_scale_factor': 10.0, 156 'static_cor_length': 100.0, 157 'lorentz_scale_factor': 3.0, 158 'dynamic_cor_length': 1.0, 159 'background': 100.0 160 }, 5.0, 100.115384615], 161 162 [{'gauss_scale_factor': 10.0, 163 'static_cor_length': 100.0, 164 'lorentz_scale_factor': 3.0, 165 'dynamic_cor_length': 1.0, 166 }, 200., 7.49981250469e-05], 167 ] -
sasmodels/models/gel_fit.py
r513efc5 r168052c 7 7 in the position of the polymer chains that ensure thermodynamic equilibrium, 8 8 and a longer distance (denoted here as $a2$ ) needed to account for the static 9 accumulations of polymer pinned down by junction points or clusters of such points.10 The latter is derived from a simple Guinier function.9 accumulations of polymer pinned down by junction points or clusters of such 10 points. The latter is derived from a simple Guinier function. 11 11 12 12 … … 40 40 --------- 41 41 42 Mitsuhiro Shibayama, Toyoichi Tanaka, Charles C Han, J. Chem. Phys. 1992, 97 (9),43 6829-684142 Mitsuhiro Shibayama, Toyoichi Tanaka, Charles C Han, 43 *J. Chem. Phys.* 1992, 97 (9), 6829-6841 44 44 45 45 Simon Mallam, Ferenc Horkay, Anne-Marie Hecht, Adrian R Rennie, Erik Geissler, 46 Macromolecules1991, 24, 543-54846 *Macromolecules* 1991, 24, 543-548 47 47 48 48 """ … … 62 62 category = "shape-independent" 63 63 64 # pylint: disable=bad-whitespace, line-too-long 64 65 # ["name", "units", default, [lower, upper], "type","description"], 65 66 parameters = [["guinier_scale", "cm^{-1}", 1.7, [-inf, inf], "", "Guinier length scale"], … … 68 69 ["fractal_exp", "", 2.0, [0, inf], "", "Fractal exponent"], 69 70 ["cor_length", "Ang", 16.0, [0, inf], "", "Correlation length"] 70 71 71 ] 72 # pylint: enable=bad-whitespace, line-too-long 72 73 73 74 source = ["gel_fit.c"] … … 92 93 'fractal_exp': 10.0, 93 94 'cor_length': 20.0 94 95 }, 0.1, 0.716532], 95 96 96 97 [{'guinier_scale': 4.0, … … 100 101 'cor_length': 20.0, 101 102 'background': 20.0, 102 }, 5.0, 20.1224653026], 103 104 ] 103 }, 5.0, 20.1224653026], 104 ] -
sasmodels/models/mass_fractal.py
r87edabf r168052c 52 52 --------- 53 53 54 D Mildner and P Hall, *J. Phys. D: Appl. Phys.*, 19 (1986) 1535-1545 Equation(9) 54 D Mildner and P Hall, *J. Phys. D: Appl. Phys.*, 55 19 (1986) 1535-1545 Equation(9) 55 56 56 57 … … 79 80 category = "shape-independent" 80 81 82 # pylint: disable=bad-whitespace, line-too-long 81 83 # ["name", "units", default, [lower, upper], "type","description"], 82 84 parameters = [["radius", "Ang", 10.0, [0.0, inf], "", "Particle radius"], 83 85 ["mass_dim", "", 1.9, [1.0, 6.0], "", "Mass fractal dimension"], 84 86 ["cutoff_length", "Ang", 100.0, [0.0, inf], "", "Cut-off length"], 85 86 87 ] 88 # pylint: enable=bad-whitespace, line-too-long 87 89 88 90 source = ["lib/sph_j1c.c", "lib/lanczos_gamma.c", "mass_fractal.c"] … … 99 101 100 102 tests = [ 101 # Accuracy tests based on content in test/utest_other_models.py102 [{'radius': 10.0,103 'mass_dim': 1.9,104 'cutoff_length': 100.0,105 }, 0.05, 279.59322],106 103 107 # Additional tests with larger range of parameters108 [{'radius': 2.0,109 'mass_dim': 3.3,110 'cutoff_length': 1.0,111 }, 0.5, 1.29016774904],104 # Accuracy tests based on content in test/utest_other_models.py 105 [{'radius': 10.0, 106 'mass_dim': 1.9, 107 'cutoff_length': 100.0, 108 }, 0.05, 279.59322], 112 109 113 [{'radius': 1.0,114 'mass_dim': 1.3,115 'cutoff_length': 1.0,116 'background': 0.8,117 }, 0.001, 1.69747015932],110 # Additional tests with larger range of parameters 111 [{'radius': 2.0, 112 'mass_dim': 3.3, 113 'cutoff_length': 1.0, 114 }, 0.5, 1.29016774904], 118 115 119 [{'radius': 1.0, 120 'mass_dim': 2.3, 121 'cutoff_length': 1.0, 122 'scale': 10.0, 123 }, 0.051, 11.6227966145], 124 ] 116 [{'radius': 1.0, 117 'mass_dim': 1.3, 118 'cutoff_length': 1.0, 119 'background': 0.8, 120 }, 0.001, 1.69747015932], 121 122 [{'radius': 1.0, 123 'mass_dim': 2.3, 124 'cutoff_length': 1.0, 125 'scale': 10.0, 126 }, 0.051, 11.6227966145], 127 ] -
sasmodels/models/mass_surface_fractal.py
r07a6700 r168052c 54 54 P Schmidt, *J Appl. Cryst.*, 24 (1991) 414-435 Equation(19) 55 55 56 A J Hurd, D W Schaefer, J E Martin, *Phys. Rev. A*, 35 (1987) 2361-2364 Equation(2) 56 A J Hurd, D W Schaefer, J E Martin, *Phys. Rev. A*, 57 35 (1987) 2361-2364 Equation(2) 57 58 58 59 """ … … 79 80 category = "shape-independent" 80 81 82 # pylint: disable=bad-whitespace, line-too-long 81 83 # ["name", "units", default, [lower, upper], "type","description"], 82 84 parameters = [["mass_dim", "", 1.8, [1e-16, 6.0], "", … … 88 90 ["primary_rg", "Ang", 4000., [0.0, inf], "", 89 91 "Primary particle radius of gyration"], 90 91 92 ] 93 # pylint: enable=bad-whitespace, line-too-long 92 94 93 95 source = ["mass_surface_fractal.c"] … … 107 109 108 110 tests = [ 109 # Accuracy tests based on content in test/utest_other_models.py110 [{'mass_dim': 1.8,111 'surface_dim': 2.3,112 'cluster_rg': 86.7,113 'primary_rg': 4000.0,114 }, 0.05, 1.77537e-05],115 111 116 # Additional tests with larger range of parameters117 [{'mass_dim': 3.3,118 'surface_dim': 1.0,119 'cluster_rg': 90.0,120 121 }, 0.001, 0.18462699016],112 # Accuracy tests based on content in test/utest_other_models.py 113 [{'mass_dim': 1.8, 114 'surface_dim': 2.3, 115 'cluster_rg': 86.7, 116 'primary_rg': 4000.0, 117 }, 0.05, 1.77537e-05], 122 118 123 [{'mass_dim': 1.3,124 'surface_dim': 1.0,125 'cluster_rg': 90.0,126 'primary_rg': 2000.0,127 'background': 0.8,128 }, 0.001, 1.16539753641],119 # Additional tests with larger range of parameters 120 [{'mass_dim': 3.3, 121 'surface_dim': 1.0, 122 'cluster_rg': 90.0, 123 'primary_rg': 4000.0, 124 }, 0.001, 0.18462699016], 129 125 130 [{'mass_dim': 2.3, 131 'surface_dim': 1.0, 132 'cluster_rg': 90.0, 133 'primary_rg': 1000.0, 134 'scale': 10.0, 135 }, 0.051, 0.000169548800377], 136 ] 126 [{'mass_dim': 1.3, 127 'surface_dim': 1.0, 128 'cluster_rg': 90.0, 129 'primary_rg': 2000.0, 130 'background': 0.8, 131 }, 0.001, 1.16539753641], 132 133 [{'mass_dim': 2.3, 134 'surface_dim': 1.0, 135 'cluster_rg': 90.0, 136 'primary_rg': 1000.0, 137 'scale': 10.0, 138 }, 0.051, 0.000169548800377], 139 ] -
sasmodels/models/polymer_excl_volume.py
r07a6700 r168052c 113 113 category = "shape-independent" 114 114 115 # pylint: disable=bad-whitespace, line-too-long 115 116 # ["name", "units", default, [lower, upper], "type", "description"], 116 117 parameters = [["rg", "Ang", 60.0, [0, inf], "", "Radius of Gyration"], 117 118 ["porod_exp", "", 3.0, [-inf, inf], "", "Porod exponent"], 118 ] 119 ] 120 # pylint: enable=bad-whitespace, line-too-long 119 121 120 122 121 def Iq(q, rg, porod_exp): 122 123 def Iq(q, 124 rg=60.0, 125 porod_exp=3.0): 123 126 """ 124 127 :param q: Input q-value (float or [float, float]) … … 131 134 o2nu = 1.0/(2.0*nu) 132 135 133 intensity = ((1.0/(nu*power(u, o2nu))) * (gamma(o2nu)*gammainc(o2nu, u) - 136 intensity = ((1.0/(nu*power(u, o2nu))) * 137 (gamma(o2nu)*gammainc(o2nu, u) - 134 138 1.0/power(u, o2nu) * gamma(porod_exp) * 135 139 gammainc(porod_exp, u))) * (q > 0) + 1.0*(q <= 0) … … 141 145 142 146 def Iqxy(qx, qy, *args): 143 iq = Iq(sqrt(qx**2 + qy**2), *args) 147 """ 148 :param qx: Input q_x-value 149 :param qy: Input q_y-value 150 :param args: Remaining arguments 151 :return: 2D-Intensity 152 """ 144 153 145 return iq154 return Iq(sqrt(qx**2 + qy**2), *args) 146 155 147 156 Iqxy.vectorized = True # Iqxy accepts an array of qx, qy values … … 158 167 159 168 tests = [ 160 161 162 163 169 # Accuracy tests based on content in test/polyexclvol_default_igor.txt 170 [{'rg': 60, 'porod_exp': 3.0}, 0.001, 0.998801], 171 [{'rg': 60, 'porod_exp': 3.0}, 0.105363, 0.0162751], 172 [{'rg': 60, 'porod_exp': 3.0}, 0.665075, 6.56261e-05], 164 173 165 166 167 168 169 170 174 # Additional tests with larger range of parameters 175 [{'rg': 10, 'porod_exp': 4.0}, 0.1, 0.723436675809], 176 [{'rg': 2.2, 'porod_exp': 22.0, 'background': 100.0}, 5.0, 100.0], 177 [{'rg': 1.1, 'porod_exp': 1, 'background': 10.0, 'scale': 1.25}, 178 20000., 10.0000712097] 179 ] -
sasmodels/models/star_polymer.py
r55b283e8 r168052c 55 55 """ 56 56 category = "shape-independent" 57 57 # pylint: disable=bad-whitespace, line-too-long 58 58 # ["name", "units", default, [lower, upper], "type","description"], 59 59 parameters = [["radius2", "Ang", 100.0, [0.0, inf], "", "Ensemble radius of gyration squared of an arm"], 60 60 ["arms", "", 3, [1.0, 6.0], "", "Number of arms in the model"], 61 62 61 ] 62 # pylint: enable=bad-whitespace, line-too-long 63 63 64 64 source = ["star_polymer.c"] … … 75 75 tests = [[{'radius2': 2.0, 76 76 'arms': 3.3, 77 77 }, 0.5, 0.850646091108], 78 78 79 79 [{'radius2': 1.0, 80 80 'arms': 2.0, 81 81 'background': 1.8, 82 83 82 }, 1.0, 2.53575888234], 83 ] -
sasmodels/models/surface_fractal.py
r07a6700 r168052c 80 80 category = "shape-independent" 81 81 82 # pylint: disable=bad-whitespace, line-too-long 82 83 # ["name", "units", default, [lower, upper], "type","description"], 83 84 parameters = [["radius", "Ang", 10.0, [0, inf], "", … … 87 88 ["cutoff_length", "Ang", 500., [0.0, inf], "", 88 89 "Cut-off Length"], 89 90 90 ] 91 # pylint: enable=bad-whitespace, line-too-long 91 92 92 93 source = ["lib/sph_j1c.c", "lib/lanczos_gamma.c", "surface_fractal.c"] … … 101 102 102 103 tests = [ 103 # Accuracy tests based on content in test/utest_other_models.py 104 [{'radius': 10.0, 'surface_dim': 2.0, 'cutoff_length': 500.0, 105 }, 0.05, 301428.65916], 104 # Accuracy tests based on content in test/utest_other_models.py 105 [{'radius': 10.0, 106 'surface_dim': 2.0, 107 'cutoff_length': 500.0, 108 }, 0.05, 301428.65916], 106 109 107 # Additional tests with larger range of parameters 108 [{'radius': 1.0, 'surface_dim': 1.0, 'cutoff_length': 10.0, 109 }, 0.332070182643, 1125.00321004], 110 # Additional tests with larger range of parameters 111 [{'radius': 1.0, 112 'surface_dim': 1.0, 113 'cutoff_length': 10.0, 114 }, 0.332070182643, 1125.00321004], 110 115 111 [{'radius': 3.5, 'surface_dim': 0.1, 'cutoff_length': 30.0, 112 'background': 0.01, 113 }, 5.0, 0.00999998891322], 116 [{'radius': 3.5, 117 'surface_dim': 0.1, 118 'cutoff_length': 30.0, 119 'background': 0.01, 120 }, 5.0, 0.00999998891322], 114 121 115 [{'radius': 3.0, 'surface_dim': 1.0, 'cutoff_length': 33.0, 116 'scale': 0.1, 117 }, 0.51, 2.50020147004], 118 ] 122 [{'radius': 3.0, 123 'surface_dim': 1.0, 124 'cutoff_length': 33.0, 125 'scale': 0.1, 126 }, 0.51, 2.50020147004], 127 ] -
sasmodels/models/two_lorentzian.py
r07a6700 r168052c 54 54 category = "shape-independent" 55 55 56 # pylint: disable=bad-whitespace, line-too-long 56 57 # ["name", "units", default, [lower, upper], "type", "description"], 57 parameters = [["lorentz_scale_1", "", 10.0, [-inf, inf], "", 58 "First power law scale factor"], 59 ["lorentz_length_1", "Ang", 100.0, [-inf, inf], "", 60 "First Lorentzian screening length"], 61 ["lorentz_exp_1", "", 3.0, [-inf, inf], "", 62 "First exponent of power law"], 63 ["lorentz_scale_2", "", 1.0, [-inf, inf], "", 64 "Second scale factor for broad Lorentzian peak"], 65 ["lorentz_length_2", "Ang", 10.0, [-inf, inf], "", 66 "Second Lorentzian screening length"], 67 ["lorentz_exp_2", "", 2.0, [-inf, inf], "", 68 "Second exponent of power law"], 69 ] 58 parameters = [["lorentz_scale_1", "", 10.0, [-inf, inf], "", "First power law scale factor"], 59 ["lorentz_length_1", "Ang", 100.0, [-inf, inf], "", "First Lorentzian screening length"], 60 ["lorentz_exp_1", "", 3.0, [-inf, inf], "", "First exponent of power law"], 61 ["lorentz_scale_2", "", 1.0, [-inf, inf], "", "Second scale factor for broad Lorentzian peak"], 62 ["lorentz_length_2", "Ang", 10.0, [-inf, inf], "", "Second Lorentzian screening length"], 63 ["lorentz_exp_2", "", 2.0, [-inf, inf], "", "Second exponent of power law"], 64 ] 65 # pylint: enable=bad-whitespace, line-too-long 70 66 71 67 72 68 def Iq(q, 73 lorentz_scale_1 ,74 lorentz_length_1 ,75 lorentz_exp_1 ,76 lorentz_scale_2 ,77 lorentz_length_2 ,78 lorentz_exp_2 ):69 lorentz_scale_1=10.0, 70 lorentz_length_1=100.0, 71 lorentz_exp_1=3.0, 72 lorentz_scale_2=1.0, 73 lorentz_length_2=10.0, 74 lorentz_exp_2=2.0): 79 75 80 76 """ … … 88 84 :return: Calculated intensity 89 85 """ 90 86 # pylint: disable=bad-whitespace 91 87 intensity = lorentz_scale_1/(1.0 + 92 88 power(q*lorentz_length_1, lorentz_exp_1)) 93 89 intensity += lorentz_scale_2/(1.0 + 94 90 power(q*lorentz_length_2, lorentz_exp_2)) 91 # pylint: enable=bad-whitespace 95 92 96 93 return intensity … … 100 97 101 98 def Iqxy(qx, qy, *args): 102 iq = Iq(sqrt(qx**2 + qy**2), *args) 99 """ 100 :param qx: Input q_x-value 101 :param qy: Input q_y-value 102 :param args: Remaining arguments 103 :return: 2D-Intensity 104 """ 103 105 104 return iq106 return Iq(sqrt(qx**2 + qy**2), *args) 105 107 106 108 Iqxy.vectorized = True # Iqxy accepts an array of qx, qy values … … 108 110 109 111 demo = dict(scale=1, background=0.1, 110 lorentz_scale_1=10, lorentz_length_1=100.0, lorentz_exp_1=3.0, 111 lorentz_scale_2=1, lorentz_length_2=10, lorentz_exp_2=2.0) 112 lorentz_scale_1=10, 113 lorentz_length_1=100.0, 114 lorentz_exp_1=3.0, 115 lorentz_scale_2=1, 116 lorentz_length_2=10, 117 lorentz_exp_2=2.0) 112 118 113 119 oldname = "TwoLorentzianModel" 114 120 oldpars = dict(background='background', 115 lorentz_scale_1='scale_1', lorentz_scale_2='scale_2', 116 lorentz_length_1='length_1', lorentz_length_2='length_2', 117 lorentz_exp_1='exponent_1', lorentz_exp_2='exponent_2') 121 lorentz_scale_1='scale_1', 122 lorentz_scale_2='scale_2', 123 lorentz_length_1='length_1', 124 lorentz_length_2='length_2', 125 lorentz_exp_1='exponent_1', 126 lorentz_exp_2='exponent_2') 118 127 119 128 tests = [ 120 # Accuracy tests based on content in test/utest_extra_models.py121 [{'lorentz_scale_1': 10.0,122 'lorentz_length_1': 100.0,123 'lorentz_exp_1': 3.0,124 'lorentz_scale_2': 1.0,125 'lorentz_length_2': 10.0,126 'lorentz_exp_2': 2.0,127 'background': 0.1,128 }, 0.001, 11.08991],129 129 130 [{'lorentz_scale_1': 10.0, 131 'lorentz_length_1': 100.0, 132 'lorentz_exp_1': 3.0, 133 'lorentz_scale_2': 1.0, 134 'lorentz_length_2': 10.0, 135 'lorentz_exp_2': 2.0, 136 'background': 0.1, 137 }, 0.150141, 0.410245], 130 # Accuracy tests based on content in test/utest_extra_models.py 131 [{'lorentz_scale_1': 10.0, 132 'lorentz_length_1': 100.0, 133 'lorentz_exp_1': 3.0, 134 'lorentz_scale_2': 1.0, 135 'lorentz_length_2': 10.0, 136 'lorentz_exp_2': 2.0, 137 'background': 0.1, 138 }, 0.001, 11.08991], 138 139 139 140 141 142 143 144 145 146 }, 0.442528, 0.148699],140 [{'lorentz_scale_1': 10.0, 141 'lorentz_length_1': 100.0, 142 'lorentz_exp_1': 3.0, 143 'lorentz_scale_2': 1.0, 144 'lorentz_length_2': 10.0, 145 'lorentz_exp_2': 2.0, 146 'background': 0.1, 147 }, 0.150141, 0.410245], 147 148 148 # Additional tests with larger range of parameters149 [{'lorentz_scale_1': 10.0,150 'lorentz_length_1': 100.0,151 'lorentz_exp_1': 3.0,152 'lorentz_scale_2': 1.0,153 'lorentz_length_2': 10.0,154 'lorentz_exp_2': 2.0,155 }, 0.000332070182643, 10.9996228107],149 [{'lorentz_scale_1': 10.0, 150 'lorentz_length_1': 100.0, 151 'lorentz_exp_1': 3.0, 152 'lorentz_scale_2': 1.0, 153 'lorentz_length_2': 10.0, 154 'lorentz_exp_2': 2.0, 155 'background': 0.1, 156 }, 0.442528, 0.148699], 156 157 157 [{'lorentz_scale_1': 0.0,158 'lorentz_length_1':0.0,159 'lorentz_exp_1':0.0,160 'lorentz_scale_2': 0.0,161 'lorentz_length_2': 0.0,162 'lorentz_exp_2':0.0,163 'background': 100.0164 }, 5.0, 100.0],158 # Additional tests with larger range of parameters 159 [{'lorentz_scale_1': 10.0, 160 'lorentz_length_1': 100.0, 161 'lorentz_exp_1': 3.0, 162 'lorentz_scale_2': 1.0, 163 'lorentz_length_2': 10.0, 164 'lorentz_exp_2': 2.0, 165 }, 0.000332070182643, 10.9996228107], 165 166 166 [{'lorentz_scale_1': 200.0, 167 'lorentz_length_1': 10.0, 168 'lorentz_exp_1': 0.1, 169 'lorentz_scale_2': 0.1, 170 'lorentz_length_2': 5.0, 171 'lorentz_exp_2': 2.0 172 }, 20000., 45.5659201896], 173 ] 167 [{'lorentz_scale_1': 0.0, 168 'lorentz_length_1': 0.0, 169 'lorentz_exp_1': 0.0, 170 'lorentz_scale_2': 0.0, 171 'lorentz_length_2': 0.0, 172 'lorentz_exp_2': 0.0, 173 'background': 100.0 174 }, 5.0, 100.0], 175 176 [{'lorentz_scale_1': 200.0, 177 'lorentz_length_1': 10.0, 178 'lorentz_exp_1': 0.1, 179 'lorentz_scale_2': 0.1, 180 'lorentz_length_2': 5.0, 181 'lorentz_exp_2': 2.0 182 }, 20000., 45.5659201896], 183 ]
Note: See TracChangeset
for help on using the changeset viewer.