1 | """ |
---|
2 | Define the resolution functions for the data. |
---|
3 | |
---|
4 | This defines classes for 1D and 2D resolution calculations. |
---|
5 | """ |
---|
6 | from scipy.special import erf |
---|
7 | from numpy import sqrt |
---|
8 | import numpy as np |
---|
9 | |
---|
10 | SLIT_SMEAR_POINTS = 500 |
---|
11 | MINIMUM_RESOLUTION = 1e-8 |
---|
12 | # TODO: Q_EXTEND_STEPS is much bigger than necessary |
---|
13 | Q_EXTEND_STEPS = 30 # number of extra q points above and below |
---|
14 | |
---|
15 | |
---|
16 | class Resolution1D(object): |
---|
17 | """ |
---|
18 | Abstract base class defining a 1D resolution function. |
---|
19 | |
---|
20 | *q* is the set of q values at which the data is measured. |
---|
21 | |
---|
22 | *q_calc* is the set of q values at which the theory needs to be evaluated. |
---|
23 | This may extend and interpolate the q values. |
---|
24 | |
---|
25 | *apply* is the method to call with I(q_calc) to compute the resolution |
---|
26 | smeared theory I(q). |
---|
27 | """ |
---|
28 | q = None |
---|
29 | q_calc = None |
---|
30 | def apply(self, Iq_calc): |
---|
31 | """ |
---|
32 | Smear *Iq_calc* by the resolution function, returning *Iq*. |
---|
33 | """ |
---|
34 | raise NotImplementedError("Subclass does not define the apply function") |
---|
35 | |
---|
36 | class Perfect1D(Resolution1D): |
---|
37 | """ |
---|
38 | Resolution function to use when there is no actual resolution smearing |
---|
39 | to be applied. It has the same interface as the other resolution |
---|
40 | functions, but returns the identity function. |
---|
41 | """ |
---|
42 | def __init__(self, q): |
---|
43 | self.q_calc = self.q = q |
---|
44 | |
---|
45 | def apply(self, Iq_calc): |
---|
46 | return Iq_calc |
---|
47 | |
---|
48 | class Pinhole1D(Resolution1D): |
---|
49 | r""" |
---|
50 | Pinhole aperture with q-dependent gaussian resolution. |
---|
51 | |
---|
52 | *q* points at which the data is measured. |
---|
53 | |
---|
54 | *dq* gaussian 1-sigma resolution at each data point. |
---|
55 | |
---|
56 | *q_calc* is the list of points to calculate, or None if this should |
---|
57 | be autogenerated from the *q,dq*. |
---|
58 | """ |
---|
59 | def __init__(self, q, q_width, q_calc=None): |
---|
60 | #TODO: maybe add min_step=np.inf |
---|
61 | #*min_step* is the minimum point spacing to use when computing the |
---|
62 | #underlying model. It should be on the order of |
---|
63 | #$\tfrac{1}{10}\tfrac{2\pi}{d_\text{max}}$ to make sure that fringes |
---|
64 | #are computed with sufficient density to avoid aliasing effects. |
---|
65 | |
---|
66 | # Protect against calls with q_width=0. The extend_q function will |
---|
67 | # not extend the q if q_width is 0, but q_width must be non-zero when |
---|
68 | # constructing the weight matrix to avoid division by zero errors. |
---|
69 | # In practice this should never be needed, since resolution should |
---|
70 | # default to Perfect1D if the pinhole geometry is not defined. |
---|
71 | self.q, self.q_width = q, q_width |
---|
72 | self.q_calc = pinhole_extend_q(q, q_width) if q_calc is None else q_calc |
---|
73 | self.weight_matrix = pinhole_resolution(self.q_calc, |
---|
74 | self.q, np.maximum(q_width, MINIMUM_RESOLUTION)) |
---|
75 | |
---|
76 | def apply(self, Iq_calc): |
---|
77 | return apply_resolution_matrix(self.weight_matrix, Iq_calc) |
---|
78 | |
---|
79 | class Slit1D(Resolution1D): |
---|
80 | """ |
---|
81 | Slit aperture with a complicated resolution function. |
---|
82 | |
---|
83 | *q* points at which the data is measured. |
---|
84 | |
---|
85 | *qx_width* slit width |
---|
86 | |
---|
87 | *qy_height* slit height |
---|
88 | """ |
---|
89 | def __init__(self, q, qx_width, qy_width): |
---|
90 | if np.isscalar(qx_width): |
---|
91 | qx_width = qx_width*np.ones(len(q)) |
---|
92 | if np.isscalar(qy_width): |
---|
93 | qy_width = qy_width*np.ones(len(q)) |
---|
94 | self.q, self.qx_width, self.qy_width = [ |
---|
95 | v.flatten() for v in (q, qx_width, qy_width)] |
---|
96 | self.q_calc = slit_extend_q(q, qx_width, qy_width) |
---|
97 | self.weight_matrix = \ |
---|
98 | slit_resolution(self.q_calc, self.q, self.qx_width, self.qy_width) |
---|
99 | |
---|
100 | def apply(self, Iq_calc): |
---|
101 | return apply_resolution_matrix(self.weight_matrix, Iq_calc) |
---|
102 | |
---|
103 | |
---|
104 | def apply_resolution_matrix(weight_matrix, Iq_calc): |
---|
105 | """ |
---|
106 | Apply the resolution weight matrix to the computed theory function. |
---|
107 | """ |
---|
108 | #print "apply shapes",Iq_calc.shape, self.weight_matrix.shape |
---|
109 | Iq = np.dot(Iq_calc[None,:], weight_matrix) |
---|
110 | #print "result shape",Iq.shape |
---|
111 | return Iq.flatten() |
---|
112 | |
---|
113 | def pinhole_resolution(q_calc, q, q_width): |
---|
114 | """ |
---|
115 | Compute the convolution matrix *W* for pinhole resolution 1-D data. |
---|
116 | |
---|
117 | Each row *W[i]* determines the normalized weight that the corresponding |
---|
118 | points *q_calc* contribute to the resolution smeared point *q[i]*. Given |
---|
119 | *W*, the resolution smearing can be computed using *dot(W,q)*. |
---|
120 | |
---|
121 | *q_calc* must be increasing. |
---|
122 | """ |
---|
123 | edges = bin_edges(q_calc) |
---|
124 | edges[edges<0.0] = 0.0 # clip edges below zero |
---|
125 | index = (q_width>0.0) # treat perfect resolution differently |
---|
126 | G = erf( (edges[:,None] - q[None,:]) / (sqrt(2.0)*q_width)[None,:] ) |
---|
127 | weights = G[1:] - G[:-1] |
---|
128 | # w = np.sum(weights, axis=0); print "w",w.shape, w |
---|
129 | weights /= np.sum(weights, axis=0)[None,:] |
---|
130 | return weights |
---|
131 | |
---|
132 | |
---|
133 | def pinhole_extend_q(q, q_width): |
---|
134 | """ |
---|
135 | Given *q* and *q_width*, find a set of sampling points *q_calc* so |
---|
136 | that each point I(q) has sufficient support from the underlying |
---|
137 | function. |
---|
138 | """ |
---|
139 | # If using min_step, you could do something like: |
---|
140 | # q_calc = np.arange(min, max, min_step) |
---|
141 | # index = np.searchsorted(q_calc, q) |
---|
142 | # q_calc[index] = q |
---|
143 | # A refinement would be to assign q to the nearest neighbour. A further |
---|
144 | # refinement would be to guard multiple q points between q_calc points, |
---|
145 | # either by removing duplicates or by only moving the values nearest the |
---|
146 | # edges. For really sparse samples, you will want to remove all remaining |
---|
147 | # points that are not within three standard deviations of a measured q. |
---|
148 | q_min, q_max = np.min(q), np.max(q) |
---|
149 | extended_min, extended_max = np.min(q - 3*q_width), np.max(q + 3*q_width) |
---|
150 | nbins_low, nbins_high = Q_EXTEND_STEPS, Q_EXTEND_STEPS |
---|
151 | if extended_min < q_min: |
---|
152 | q_low = np.linspace(extended_min, q_min, nbins_low+1)[:-1] |
---|
153 | q_low = q_low[q_low>0.0] |
---|
154 | else: |
---|
155 | q_low = [] |
---|
156 | if extended_max > q_max: |
---|
157 | q_high = np.linspace(q_max, extended_max, nbins_high+1)[:-1] |
---|
158 | else: |
---|
159 | q_high = [] |
---|
160 | return np.concatenate([q_low, q, q_high]) |
---|
161 | |
---|
162 | |
---|
163 | def slit_resolution(q_calc, q, qx_width, qy_width): |
---|
164 | edges = bin_edges(q_calc) # Note: requires q > 0 |
---|
165 | edges[edges<0.0] = 0.0 # clip edges below zero |
---|
166 | qy_min, qy_max = 0.0, edges[-1] |
---|
167 | |
---|
168 | # Make q_calc into a row vector, and q, qx_width, qy_width into columns |
---|
169 | # Make weights match [ q_calc X q ] |
---|
170 | weights = np.zeros((len(q),len(q_calc)),'d') |
---|
171 | q_calc = q_calc[None,:] |
---|
172 | q, qx_width, qy_width, edges = [ |
---|
173 | v[:,None] for v in (q, qx_width, qy_width, edges)] |
---|
174 | |
---|
175 | # Loop for width (height is analytical). |
---|
176 | # Condition: height >>> width, otherwise, below is not accurate enough. |
---|
177 | # Smear weight numerical iteration for width>0 when height>0. |
---|
178 | # When width = 0, the numerical iteration will be skipped. |
---|
179 | # The resolution calculation for the height is done by direct integration, |
---|
180 | # assuming the I(q'=sqrt(q_j^2-(q+shift_w)^2)) is constant within |
---|
181 | # a q' bin, [q_high, q_low]. |
---|
182 | # In general, this weight numerical iteration for width>0 might be a rough |
---|
183 | # approximation, but it must be good enough when height >>> width. |
---|
184 | E_sq = edges**2 |
---|
185 | y_points = SLIT_SMEAR_POINTS if np.any(qy_width>0) else 1 |
---|
186 | qy_step = 0 if y_points == 1 else qy_width/(y_points-1) |
---|
187 | for k in range(-y_points+1,y_points): |
---|
188 | qy = np.clip(q + qy_step*k, qy_min, qy_max) |
---|
189 | qx_low = qy |
---|
190 | qx_high = sqrt(qx_low**2 + qx_width**2) |
---|
191 | in_x = (q_calc>=qx_low)*(q_calc<=qx_high) |
---|
192 | qy_sq = qy**2 |
---|
193 | weights += (sqrt(E_sq[1:]-qy_sq) - sqrt(qy_sq - E_sq[:-1]))*in_x |
---|
194 | w = np.sum(weights, axis=1); print "w",w.shape, w |
---|
195 | weights /= np.sum(weights, axis=1)[:,None] |
---|
196 | return weights |
---|
197 | |
---|
198 | |
---|
199 | def slit_extend_q(q, qx_width, qy_width): |
---|
200 | return q |
---|
201 | |
---|
202 | |
---|
203 | def bin_edges(x): |
---|
204 | """ |
---|
205 | Determine bin edges from bin centers, assuming that edges are centered |
---|
206 | between the bins. |
---|
207 | |
---|
208 | Note: this uses the arithmetic mean, which may not be appropriate for |
---|
209 | log-scaled data. |
---|
210 | """ |
---|
211 | if len(x) < 2 or (np.diff(x)<0).any(): |
---|
212 | raise ValueError("Expected bins to be an increasing set") |
---|
213 | edges = np.hstack([ |
---|
214 | x[0] - 0.5*(x[1] - x[0]), # first point minus half first interval |
---|
215 | 0.5*(x[1:] + x[:-1]), # mid points of all central intervals |
---|
216 | x[-1] + 0.5*(x[-1] - x[-2]), # last point plus half last interval |
---|
217 | ]) |
---|
218 | return edges |
---|
219 | |
---|
220 | ############################################################################ |
---|
221 | # unit tests |
---|
222 | ############################################################################ |
---|
223 | import unittest |
---|
224 | |
---|
225 | class ResolutionTest(unittest.TestCase): |
---|
226 | |
---|
227 | def setUp(self): |
---|
228 | self.x = 0.001*np.arange(1,11) |
---|
229 | self.y = self.Iq(self.x) |
---|
230 | |
---|
231 | def Iq(self, q): |
---|
232 | "Linear function for resolution unit test" |
---|
233 | return 12.0 - 1000.0*q |
---|
234 | |
---|
235 | def test_perfect(self): |
---|
236 | """ |
---|
237 | Perfect resolution and no smearing. |
---|
238 | """ |
---|
239 | resolution = Perfect1D(self.x) |
---|
240 | Iq_calc = self.Iq(resolution.q_calc) |
---|
241 | output = resolution.apply(Iq_calc) |
---|
242 | np.testing.assert_equal(output, self.y) |
---|
243 | |
---|
244 | def test_slit_zero(self): |
---|
245 | """ |
---|
246 | Slit smearing with perfect resolution. |
---|
247 | """ |
---|
248 | resolution = Slit1D(self.x, 0., 0.) |
---|
249 | Iq_calc = self.Iq(resolution.q_calc) |
---|
250 | output = resolution.apply(Iq_calc) |
---|
251 | np.testing.assert_equal(output, self.y) |
---|
252 | |
---|
253 | @unittest.skip("slit smearing is still broken") |
---|
254 | def test_slit(self): |
---|
255 | """ |
---|
256 | Slit smearing with height 0.005 |
---|
257 | """ |
---|
258 | resolution = Slit1D(self.x, 0., 0.005) |
---|
259 | Iq_calc = self.Iq(resolution.q_calc) |
---|
260 | output = resolution.apply(Iq_calc) |
---|
261 | # The following commented line was the correct output for even bins [see smearer.cpp for details] |
---|
262 | #answer = [ 9.666, 9.056, 8.329, 7.494, 6.642, 5.721, 4.774, 3.824, 2.871, 2. ] |
---|
263 | answer = [ 9.0618, 8.6401, 8.1186, 7.1391, 6.1528, 5.5555, 4.5584, 3.5606, 2.5623, 2. ] |
---|
264 | np.testing.assert_allclose(output, answer, atol=1e-4) |
---|
265 | |
---|
266 | def test_pinhole_zero(self): |
---|
267 | """ |
---|
268 | Pinhole smearing with perfect resolution |
---|
269 | """ |
---|
270 | resolution = Pinhole1D(self.x, 0.0*self.x) |
---|
271 | Iq_calc = self.Iq(resolution.q_calc) |
---|
272 | output = resolution.apply(Iq_calc) |
---|
273 | np.testing.assert_equal(output, self.y) |
---|
274 | |
---|
275 | def test_pinhole(self): |
---|
276 | """ |
---|
277 | Pinhole smearing with dQ = 0.001 [Note: not dQ/Q = 0.001] |
---|
278 | """ |
---|
279 | resolution = Pinhole1D(self.x, 0.001*np.ones_like(self.x), |
---|
280 | q_calc=self.x) |
---|
281 | Iq_calc = 12.0-1000.0*resolution.q_calc |
---|
282 | output = resolution.apply(Iq_calc) |
---|
283 | answer = [ 10.44785079, 9.84991299, 8.98101708, |
---|
284 | 7.99906585, 6.99998311, 6.00001689, |
---|
285 | 5.00093415, 4.01898292, 3.15008701, 2.55214921] |
---|
286 | np.testing.assert_allclose(output, answer, atol=1e-8) |
---|
287 | |
---|
288 | # Q, dQ, I(Q) for Igor default sphere model. |
---|
289 | # combines CMSphere5.txt and CMSphere5smearsphere.txt from sasview/tests |
---|
290 | # TODO: move test data into its own file? |
---|
291 | SPHERE_RESOLUTION_TEST_DATA = """\ |
---|
292 | 0.001278 0.0002847 2538.41176383 |
---|
293 | 0.001562 0.0002905 2536.91820405 |
---|
294 | 0.001846 0.0002956 2535.13182479 |
---|
295 | 0.002130 0.0003017 2533.06217813 |
---|
296 | 0.002414 0.0003087 2530.70378586 |
---|
297 | 0.002698 0.0003165 2528.05024192 |
---|
298 | 0.002982 0.0003249 2525.10408349 |
---|
299 | 0.003266 0.0003340 2521.86667499 |
---|
300 | 0.003550 0.0003437 2518.33907750 |
---|
301 | 0.003834 0.0003539 2514.52246995 |
---|
302 | 0.004118 0.0003646 2510.41798319 |
---|
303 | 0.004402 0.0003757 2506.02690988 |
---|
304 | 0.004686 0.0003872 2501.35067884 |
---|
305 | 0.004970 0.0003990 2496.38678318 |
---|
306 | 0.005253 0.0004112 2491.16237596 |
---|
307 | 0.005537 0.0004237 2485.63911673 |
---|
308 | 0.005821 0.0004365 2479.83657083 |
---|
309 | 0.006105 0.0004495 2473.75676948 |
---|
310 | 0.006389 0.0004628 2467.40145990 |
---|
311 | 0.006673 0.0004762 2460.77293372 |
---|
312 | 0.006957 0.0004899 2453.86724627 |
---|
313 | 0.007241 0.0005037 2446.69623838 |
---|
314 | 0.007525 0.0005177 2439.25775219 |
---|
315 | 0.007809 0.0005318 2431.55421398 |
---|
316 | 0.008093 0.0005461 2423.58785521 |
---|
317 | 0.008377 0.0005605 2415.36158137 |
---|
318 | 0.008661 0.0005750 2406.87009473 |
---|
319 | 0.008945 0.0005896 2398.12841186 |
---|
320 | 0.009229 0.0006044 2389.13360806 |
---|
321 | 0.009513 0.0006192 2379.88958042 |
---|
322 | 0.009797 0.0006341 2370.39776774 |
---|
323 | 0.010080 0.0006491 2360.69528793 |
---|
324 | 0.010360 0.0006641 2350.85169027 |
---|
325 | 0.010650 0.0006793 2340.42023633 |
---|
326 | 0.010930 0.0006945 2330.11206013 |
---|
327 | 0.011220 0.0007097 2319.20109972 |
---|
328 | 0.011500 0.0007251 2308.43503981 |
---|
329 | 0.011780 0.0007404 2297.44820179 |
---|
330 | 0.012070 0.0007558 2285.83853677 |
---|
331 | 0.012350 0.0007713 2274.41290746 |
---|
332 | 0.012640 0.0007868 2262.36219581 |
---|
333 | 0.012920 0.0008024 2250.51169731 |
---|
334 | 0.013200 0.0008180 2238.45596231 |
---|
335 | 0.013490 0.0008336 2225.76495666 |
---|
336 | 0.013770 0.0008493 2213.29618391 |
---|
337 | 0.014060 0.0008650 2200.19110751 |
---|
338 | 0.014340 0.0008807 2187.34050325 |
---|
339 | 0.014620 0.0008965 2174.30529864 |
---|
340 | 0.014910 0.0009123 2160.61632548 |
---|
341 | 0.015190 0.0009281 2147.21038112 |
---|
342 | 0.015470 0.0009440 2133.62023580 |
---|
343 | 0.015760 0.0009598 2119.37907426 |
---|
344 | 0.016040 0.0009757 2105.45234903 |
---|
345 | 0.016330 0.0009916 2090.86319102 |
---|
346 | 0.016610 0.0010080 2076.60576032 |
---|
347 | 0.016890 0.0010240 2062.19214565 |
---|
348 | 0.017180 0.0010390 2047.10550219 |
---|
349 | 0.017460 0.0010550 2032.38715621 |
---|
350 | 0.017740 0.0010710 2017.52560123 |
---|
351 | 0.018030 0.0010880 2001.99124318 |
---|
352 | 0.018310 0.0011040 1986.84662060 |
---|
353 | 0.018600 0.0011200 1971.03389745 |
---|
354 | 0.018880 0.0011360 1955.61395119 |
---|
355 | 0.019160 0.0011520 1940.08291563 |
---|
356 | 0.019450 0.0011680 1923.87672225 |
---|
357 | 0.019730 0.0011840 1908.10656374 |
---|
358 | 0.020020 0.0012000 1891.66297192 |
---|
359 | 0.020300 0.0012160 1875.66789021 |
---|
360 | 0.020580 0.0012320 1859.56357196 |
---|
361 | 0.020870 0.0012490 1842.79468290 |
---|
362 | 0.021150 0.0012650 1826.50064489 |
---|
363 | 0.021430 0.0012810 1810.11533702 |
---|
364 | 0.021720 0.0012970 1793.06840882 |
---|
365 | 0.022000 0.0013130 1776.51153580 |
---|
366 | 0.022280 0.0013290 1759.87201249 |
---|
367 | 0.022570 0.0013460 1742.57354412 |
---|
368 | 0.022850 0.0013620 1725.79397319 |
---|
369 | 0.023140 0.0013780 1708.35831550 |
---|
370 | 0.023420 0.0013940 1691.45256069 |
---|
371 | 0.023700 0.0014110 1674.48561783 |
---|
372 | 0.023990 0.0014270 1656.86525366 |
---|
373 | 0.024270 0.0014430 1639.79847285 |
---|
374 | 0.024550 0.0014590 1622.68887088 |
---|
375 | 0.024840 0.0014760 1604.96421100 |
---|
376 | 0.025120 0.0014920 1587.85768129 |
---|
377 | 0.025410 0.0015080 1569.99297335 |
---|
378 | 0.025690 0.0015240 1552.84580279 |
---|
379 | 0.025970 0.0015410 1535.54074115 |
---|
380 | 0.026260 0.0015570 1517.75249337 |
---|
381 | 0.026540 0.0015730 1500.40115023 |
---|
382 | 0.026820 0.0015900 1483.03632237 |
---|
383 | 0.027110 0.0016060 1465.05942429 |
---|
384 | 0.027390 0.0016220 1447.67682181 |
---|
385 | 0.027670 0.0016390 1430.46495191 |
---|
386 | 0.027960 0.0016550 1412.49232282 |
---|
387 | 0.028240 0.0016710 1395.13182318 |
---|
388 | 0.028520 0.0016880 1377.93439837 |
---|
389 | 0.028810 0.0017040 1359.99528971 |
---|
390 | 0.029090 0.0017200 1342.67274512 |
---|
391 | 0.029370 0.0017370 1325.55375609 |
---|
392 | """ |
---|
393 | |
---|
394 | class Pinhole1DTest(unittest.TestCase): |
---|
395 | |
---|
396 | def setUp(self): |
---|
397 | # sample q, dq, I(q) calculated by NIST Igor SANS package |
---|
398 | self.data = np.loadtxt(SPHERE_RESOLUTION_TEST_DATA.split('\n')).T |
---|
399 | self.pars = dict(scale=1.0, background=0.01, radius=60.0, |
---|
400 | solvent_sld=6.3, sld=1) |
---|
401 | |
---|
402 | def Iq(self, q, dq): |
---|
403 | from sasmodels import core |
---|
404 | from sasmodels.models import sphere |
---|
405 | |
---|
406 | model = core.load_model(sphere) |
---|
407 | resolution = Pinhole1D(q, dq) |
---|
408 | kernel = core.make_kernel(model, [resolution.q_calc]) |
---|
409 | Iq_calc = core.call_kernel(kernel, self.pars) |
---|
410 | result = resolution.apply(Iq_calc) |
---|
411 | return result |
---|
412 | |
---|
413 | def test_sphere(self): |
---|
414 | """ |
---|
415 | Compare pinhole resolution smearing with NIST Igor SANS |
---|
416 | """ |
---|
417 | q, dq, answer = self.data |
---|
418 | output = self.Iq(q, dq) |
---|
419 | np.testing.assert_allclose(output, answer, rtol=0.006) |
---|
420 | |
---|
421 | #TODO: is all sas data sampled densely enough to support resolution calcs? |
---|
422 | @unittest.skip("suppress sparse data test; not supported by current code") |
---|
423 | def test_sphere_sparse(self): |
---|
424 | """ |
---|
425 | Compare pinhole resolution smearing with NIST Igor SANS on sparse data |
---|
426 | """ |
---|
427 | q, dq, answer = self.data[:, ::20] # Take every nth point |
---|
428 | output = self.Iq(q, dq) |
---|
429 | np.testing.assert_allclose(output, answer, rtol=0.006) |
---|
430 | |
---|
431 | |
---|
432 | def main(): |
---|
433 | """ |
---|
434 | Run tests given is sys.argv. |
---|
435 | |
---|
436 | Returns 0 if success or 1 if any tests fail. |
---|
437 | """ |
---|
438 | import sys |
---|
439 | import xmlrunner |
---|
440 | |
---|
441 | suite = unittest.TestSuite() |
---|
442 | suite.addTest(unittest.defaultTestLoader.loadTestsFromModule(sys.modules[__name__])) |
---|
443 | |
---|
444 | runner = xmlrunner.XMLTestRunner(output='logs') |
---|
445 | result = runner.run(suite) |
---|
446 | return 1 if result.failures or result.errors else 0 |
---|
447 | |
---|
448 | |
---|
449 | ############################################################################ |
---|
450 | # usage demo |
---|
451 | ############################################################################ |
---|
452 | |
---|
453 | def _eval_demo_1d(resolution, title): |
---|
454 | from sasmodels import core |
---|
455 | from sasmodels.models import cylinder |
---|
456 | ## or alternatively: |
---|
457 | # cylinder = core.load_model_definition('cylinder') |
---|
458 | model = core.load_model(cylinder) |
---|
459 | |
---|
460 | kernel = core.make_kernel(model, [resolution.q_calc]) |
---|
461 | Iq_calc = core.call_kernel(kernel, {'length':210, 'radius':500}) |
---|
462 | Iq = resolution.apply(Iq_calc) |
---|
463 | |
---|
464 | import matplotlib.pyplot as plt |
---|
465 | plt.loglog(resolution.q_calc, Iq_calc, label='unsmeared') |
---|
466 | plt.loglog(resolution.q, Iq, label='smeared', hold=True) |
---|
467 | plt.legend() |
---|
468 | plt.title(title) |
---|
469 | plt.xlabel("Q (1/Ang)") |
---|
470 | plt.ylabel("I(Q) (1/cm)") |
---|
471 | |
---|
472 | def demo_pinhole_1d(): |
---|
473 | q = np.logspace(-3,-1,400) |
---|
474 | dq = 0.1*q |
---|
475 | resolution = Pinhole1D(q, dq) |
---|
476 | _eval_demo_1d(resolution, title="10% dQ/Q Pinhole Resolution") |
---|
477 | |
---|
478 | def demo_slit_1d(): |
---|
479 | q = np.logspace(-3,-1,400) |
---|
480 | qx_width = 0.005 |
---|
481 | qy_width = 0.0 |
---|
482 | resolution = Slit1D(q, qx_width, qy_width) |
---|
483 | _eval_demo_1d(resolution, title="0.005 Qx Slit Resolution") |
---|
484 | |
---|
485 | def demo(): |
---|
486 | import matplotlib.pyplot as plt |
---|
487 | plt.subplot(121) |
---|
488 | demo_pinhole_1d() |
---|
489 | plt.subplot(122) |
---|
490 | demo_slit_1d() |
---|
491 | plt.show() |
---|
492 | |
---|
493 | |
---|
494 | if __name__ == "__main__": |
---|
495 | #demo() |
---|
496 | main() |
---|
497 | |
---|
498 | |
---|