diffpy.structure/src/diffpy/structure/expansion/makeellipsoid.py at 5f64982363999d96543df9934fc6e13250454032 · stevenhua0320/diffpy.structure · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
#!/usr/bin/env python
##############################################################################
#
# diffpy.structure  by DANSE Diffraction group
#                   Simon J. L. Billinge
#                   (c) 2008 trustees of the Michigan State University.
#                   All rights reserved.
#
# File coded by:    Chris Farrow
#
# See AUTHORS.txt for a list of people who contributed.
# See LICENSE_DANSE.txt for license information.
#
##############################################################################
"""Make a spheroid nanoparticle from a template structure."""

from math import ceil

from numpy import array

from diffpy.structure import Structure
from diffpy.structure.expansion.shapeutils import find_center
from diffpy.utils._deprecator import build_deprecation_message, deprecated

base = "diffpy.structure"
removal_version = "4.0.0"
makeSphere_deprecation_msg = build_deprecation_message(
    base,
    "makeSphere",
    "make_sphere",
    removal_version,
)
makeEllipsoid_deprecation_msg = build_deprecation_message(
    base,
    "makeEllipsoid",
    "make_ellipsoid",
    removal_version,
)


@deprecated(makeSphere_deprecation_msg)
def makeSphere(S, radius):
    """This function has been deprecated and will be removed in version
    4.0.0.

    Please use diffpy.structure.make_sphere instead.
    """
    return make_sphere(S, radius)


def make_sphere(S, radius):
    """Create a spherical nanoparticle.

    Parameters
    ----------
    S : Structure
        A `Structure` instance.
    radius : float
        Primary equatorial radius (along x-axis).

    Returns
    -------
    Structure
        A new `Structure` instance.
    """
    return make_ellipsoid(S, radius)


@deprecated(makeEllipsoid_deprecation_msg)
def makeEllipsoid(S, a, b=None, c=None):
    """This function has been deprecated and will be removed in version
    4.0.0.

    Please use diffpy.structure.make_ellipsoid instead.
    """
    return make_ellipsoid(S, a, b, c)


def make_ellipsoid(S, a, b=None, c=None):
    """Cut a `Structure` out of another one.

    Parameters
    ----------
    S : Structure
        A `Structure` instance.
    a : float
        Primary equatorial radius (along x-axis).
    b : float, Optional
        Secondary equatorial radius (along y-axis). If `b` is ``None``
        (default), then it is set equal to `a`.
    c : float, Optional
        Polar radius (along z-axis). If `c` is ``None`` (default), then it is
        set equal to `a`.

    Returns
    -------
    Structure :
        A new `Structure` instance.
    """
    if b is None:
        b = a
    if c is None:
        c = a
    sabc = array([a, b, c])

    # Create a supercell large enough for the ellipsoid
    frac = S.lattice.fractional(sabc)
    # FIXME - this looks fishy for non-orthogonal lattices
    mno = max(ceil(2 * xi) for xi in frac) * array([1, 1, 1])
    # Make the supercell
    from diffpy.structure.expansion import supercell

    newS = supercell(S, mno)
    lat = newS.lattice

    # Find the central atom
    ncenter = find_center(newS)

    cxyz = lat.cartesian(newS[ncenter].xyz)

    delList = []
    N = len(newS)
    j = N
    for i in range(N):
        j -= 1

        # Calculate (x/a)**2 + (y/b)**2 + (z/c)**2
        xyz = lat.cartesian(newS[j].xyz)
        darray = ((xyz - cxyz) / sabc) ** 2
        d = sum(darray) ** 0.5

        # Discard atom if (x/a)**2 + (y/b)**2 + (z/c)**2 > 1
        if d > 1:
            delList.append(j)

    for i in delList:
        newS.pop(i)

    return newS


# ----------------------------------------------------------------------------

if __name__ == "__main__":
    import os.path

    datadir = "../../tests/testdata"
    S = Structure()
    S.read(os.path.join(datadir, "CdSe_bulk.stru"), "pdffit")
    newS = make_ellipsoid(S, 12)
    newS.write("CdSe_d24.stru", "pdffit")
    newS = make_ellipsoid(S, 20, 10, 10)
    newS.write("CdSe_a20_b10_c10.stru", "pdffit")
    newS = make_ellipsoid(S, 20, 15, 10)
    newS.write("CdSe_a20_b15_c10.stru", "pdffit")
    S = Structure()
    S.read(os.path.join(datadir, "Ni.stru"), "pdffit")
    newS = make_ellipsoid(S, 10)
    newS.write("Ni_d20.stru", "pdffit")
    newS = make_ellipsoid(S, 20, 4)
    newS.write("Ni_a20_b4_c20.stru", "pdffit")
    newS = make_ellipsoid(S, 20, 15, 10)
    newS.write("Ni_a20_b15_c10.stru", "pdffit")