Changeset 2ab331f in sasmodels

Timestamp:

Jan 11, 2018 11:00:36 AM (7 years ago)

Author:

Paul Kienzle <pkienzle@…>

Branches:

master, core_shell_microgels, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests

Children:

2a7e20e

Parents:

226473d

Message:

play with non-dilute cylinder models

File:

• explore/realspace.py (modified) (5 diffs)

explore/realspace.py

@@ -393 +393 @@
         pass
     n = len(points)
-    index = np.random.choice(n, size=1000) if n > 1000 else slice(None, None)
+    #print("len points", n)
+    index = np.random.choice(n, size=500) if n > 500 else slice(None, None)
     ax.scatter(points[index, 0], points[index, 1], points[index, 2], c=rho[index])
     #low, high = points.min(axis=0), points.max(axis=0)
@@ -476 +477 @@
     q = np.logspace(-3, 0, 200)
     r = shape.r_bins(q, r_step=0.01)
-    sampling_density = 15000 / shape.volume()
+    sampling_density = 50000 / shape.volume()
     rho, points = shape.sample(sampling_density)
     Pr = calc_Pr(r, rho-rho_solvent, points)
@@ -489 +490 @@
 def check_shape_2d(shape, fn=None, view=(0, 0, 0)):
     rho_solvent = 0
-    qx = qy = np.linspace(-1, 1, 100)
+    nq, qmax = 100, 0.1
+    nq, qmax = 100, 0.03
+    qx = np.linspace(0.0, qmax, nq)
+    qy = np.linspace(0.0, qmax, nq)
     Qx, Qy = np.meshgrid(qx, qy)
-    sampling_density = 50000 / shape.volume()
+    sampling_density = 15000 / shape.volume()
+    #t0 = time.time()
     rho, points = shape.sample(sampling_density)
+    #print("sample time", time.time() - t0)
+    t0 = time.time()
     Iqxy = calc_Iqxy(Qx, Qy, rho, points, view=view)
+    print("calc time", time.time() - t0)
+    theory = fn(Qx, Qy) if fn is not None else None
     Iqxy += 0.001 * Iqxy.max()
-    theory = fn(Qx, Qy)+0.001 if fn is not None else None
+    if theory is not None:
+        theory += 0.001 * theory.max()
 
     import pylab
+    #plot_points(rho, points); pylab.figure()
     plot_calc_2d(qx, qy, Iqxy, theory=theory)
     pylab.show()
@@ -508 +519 @@
 def check_cylinder_2d(radius=25, length=125, rho=2., view=(0, 0, 0)):
     shape = EllipticalCylinder(radius, radius, length, rho)
+    fn = lambda qx, qy, view=view: cylinder_Iqxy(qx, qy, radius, length, view=view)
+    check_shape_2d(shape, fn, view=view)
+
+def check_cylinder_2d_lattice(radius=25, length=125, rho=2.,
+                              view=(0, 0, 0)):
+    nx, dx = 1, 2*radius
+    ny, dy = 30, 2*radius
+    nz, dz = 30, length
+    dx, dy, dz = 2*dx, 2*dy, 2*dz
+    def center(*args):
+        sigma = 0.333
+        space = 2
+        return [(space*n+np.random.randn()*sigma)*x for n, x in args]
+    shapes = [EllipticalCylinder(radius, radius, length, rho,
+                                 #center=(ix*dx, iy*dy, iz*dz)
+                                 orientation=np.random.randn(3)*0,
+                                 center=center((ix, dx), (iy, dy), (iz, dz))
+                                )
+              for ix in range(nx)
+              for iy in range(ny)
+              for iz in range(nz)]
+    shape = Composite(shapes)
     fn = lambda qx, qy, view=view: cylinder_Iqxy(qx, qy, radius, length, view=view)
     check_shape_2d(shape, fn, view=view)
@@ -531 +564 @@
     check_cylinder(radius=10, length=40)
     #check_cylinder_2d(radius=10, length=40, view=(90,30,0))
+    #check_cylinder_2d_lattice(radius=10, length=50, view=(90,30,0))
     #check_sphere()
     #check_csbox()