- Timestamp:
- Sep 12, 2018 5:13:27 PM (6 years ago)
- Branches:
- master
- Children:
- 01dba26
- Parents:
- d0fdba2
- Location:
- example/weights
- Files:
-
- 3 added
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
example/weights/cyclic_gaussian.py
r36873d1 r55e82f0 7 7 r""" 8 8 Cyclic gaussian dispersion on orientation. 9 9 10 10 .. math: 11 11 … … 13 13 14 14 This provides a close match to the gaussian distribution for 15 low angles (with $\sin \theta \approx \theta$), but the tails 16 are limited to $\pm 90^\circ$. For $\sigma$ large the 17 distribution is approximately uniform. The usual polar coordinate 15 low angles, but the tails are limited to $\pm 90^\circ$. For $\sigma$ 16 large the distribution is approximately uniform. The usual polar coordinate 18 17 projection applies, with $\theta$ weights scaled by $\cos \theta$ 19 18 and $\phi$ weights unscaled. 20 19 21 This is closely related to a Maier-Saupe distribution with order 22 parameter $P_2$ and appropriate scaling constants, and changes 23 between $\sin$ and $\cos$ as appropriate for the coordinate system 24 representation. 20 This is eqivalent to a Maier-Saupe distribution with order 21 parameter $P_2 = 1/(2 \sigma^2)$, with $\sigma$ in radians. 25 22 """ 26 23 type = "cyclic_gaussian" … … 29 26 # Note: center is always zero for orientation distributions 30 27 def _weights(self, center, sigma, lb, ub): 31 # Convert sigma in degrees to the approximately equivalent Maier-Saupe "a"28 # Convert sigma in degrees to radians 32 29 sigma = radians(sigma) 33 a = -0.5/sigma**234 30 35 # Limit width to +/-90 degrees; use an open interval since the 36 # pattern at +90 is the same as that at -90. 31 # Limit width to +/- 90 degrees 37 32 width = min(self.nsigmas*sigma, pi/2) 38 x = np.linspace(-width, width, self.npts+2)[1:-1] 33 x = np.linspace(-width, width, self.npts) 34 35 # Truncate the distribution in case the parameter value is limited 36 x[(x >= radians(lb)) & (x <= radians(ub))] 39 37 40 38 # Return orientation in degrees with Maier-Saupe weights 41 return degrees(x), exp(a*sin(x)**2) 42 39 return degrees(x), exp(-0.5*sin(x)**2/sigma**2)
Note: See TracChangeset
for help on using the changeset viewer.