From b6e2f7875029d1a5a0f6f7a5761b5f3af1f7b448 Mon Sep 17 00:00:00 2001 From: Arun Isaac Date: Fri, 26 Mar 2021 14:43:07 +0530 Subject: Simplify directory structure. * setup.cfg: Specify the sambal package explicitly. * src/sambal: Move to sambal. --- sambal/sambal.py | 68 ++++++++++++++++++++++++++++++++++++++++++++++++++ setup.cfg | 7 +----- src/sambal/__init__.py | 0 src/sambal/sambal.py | 65 ----------------------------------------------- 4 files changed, 69 insertions(+), 71 deletions(-) create mode 100644 sambal/sambal.py delete mode 100644 src/sambal/__init__.py delete mode 100644 src/sambal/sambal.py diff --git a/sambal/sambal.py b/sambal/sambal.py new file mode 100644 index 0000000..f4d64f4 --- /dev/null +++ b/sambal/sambal.py @@ -0,0 +1,68 @@ +# sambal --- Sample balls, spheres, spherical caps +# Copyright © 2021 Arun I +# Copyright © 2021 Murugesan Venkatapathi +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, but +# WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +# General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see +# . + +from numpy import cos, dot, empty, log, sin, sqrt, pi +from numpy.random import randn, random +from numpy.linalg import norm + +def random_on_sphere(dim): + """Return a random vector uniformly distributed on the unit sphere.""" + x = randn(dim) + return x / norm(x) + +def rotate_from_nth_canonical(x, axis): + """Rotate vector from around the nth canonical axis to the given axis. + + """ + xn = x[-1] + axisn = axis[-1] + if axisn != 1: + b = norm(axis[:-1]) + a = (dot(x, axis) - xn*axisn) / b + s = sqrt(1 - axisn**2) + x = x + (xn*s + a*(axisn - 1))/b * axis + x[-1] = x[-1] + xn*(axisn - 1) - a*s \ + - axisn*(xn*s + a*(axisn - 1))/b + return x + +def random_on_disk(axis, planar_angle): + """Return a random vector uniformly distributed on the periphery of a +disk. + + """ + dim = axis.size + x = empty(dim) + x[:-1] = sin(planar_angle) * random_on_sphere(dim - 1) + x[-1] = cos(planar_angle) + return rotate_from_nth_canonical(x, axis) + +def random_on_cap(axis, maximum_planar_angle): + """Return a random vector uniformly distributed on a spherical +cap. The random planar angle is generated using rejection sampling. + + """ + # We apply the log function just to prevent the floats from + # underflowing. + dim = axis.size + box_height = (dim-2)*log(sin(min(maximum_planar_angle, pi/2))) + while True: + theta = maximum_planar_angle*random() + f = box_height + log(random()) + if f < (dim-2)*log(sin(theta)): + break + return random_on_disk(axis, theta) diff --git a/setup.cfg b/setup.cfg index 604dfa7..984bd38 100644 --- a/setup.cfg +++ b/setup.cfg @@ -13,9 +13,4 @@ classifiers = Operating System :: OS Independent [options] -package_dir = - = src -packages = find: - -[options.packages.find] -where = src \ No newline at end of file +packages = sambal diff --git a/src/sambal/__init__.py b/src/sambal/__init__.py deleted file mode 100644 index e69de29..0000000 diff --git a/src/sambal/sambal.py b/src/sambal/sambal.py deleted file mode 100644 index c42ea44..0000000 --- a/src/sambal/sambal.py +++ /dev/null @@ -1,65 +0,0 @@ -# sambal --- Sample balls, spheres, spherical caps -# Copyright © 2021 Arun I -# Copyright © 2021 Murugesan Venkatapathi -# -# This program is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 3 of the License, or -# (at your option) any later version. -# -# This program is distributed in the hope that it will be useful, but -# WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -# General Public License for more details. -# -# You should have received a copy of the GNU General Public License -# along with this program. If not, see -# . - -from numpy import cos, dot, empty, log, sin, sqrt, pi -from numpy.random import randn, random -from numpy.linalg import norm - -def random_on_sphere(dim): - """Return a random vector uniformly distributed on the unit sphere.""" - x = randn(dim) - return x / norm(x) - -def rotate_from_nth_canonical(x, axis): - """Rotate vector from around the nth canonical axis to the given axis. - """ - xn = x[-1] - axisn = axis[-1] - if axisn != 1: - b = norm(axis[:-1]) - a = (dot(x, axis) - xn*axisn) / b - s = sqrt(1 - axisn**2) - x = x + (xn*s + a*(axisn - 1))/b * axis - x[-1] = x[-1] + xn*(axisn - 1) - a*s \ - - axisn*(xn*s + a*(axisn - 1))/b - return x - -def random_on_disk(axis, planar_angle): - """Return a random vector uniformly distributed on the periphery of a -disk.""" - dim = axis.size - x = empty(dim) - x[:-1] = sin(planar_angle) * random_on_sphere(dim - 1) - x[-1] = cos(planar_angle) - return rotate_from_nth_canonical(x, axis) - -def random_on_cap(axis, maximum_planar_angle): - """Return a random vector uniformly distributed on a spherical -cap. The random planar angle is generated using rejection sampling. - - """ - # We apply the log function just to prevent the floats from - # underflowing. - dim = axis.size - box_height = (dim-2)*log(sin(min(maximum_planar_angle, pi/2))) - while True: - theta = maximum_planar_angle*random() - f = box_height + log(random()) - if f < (dim-2)*log(sin(theta)): - break - return random_on_disk(axis, theta) -- cgit v1.2.3