Routines for finding vortex in Bose-Einstein condensate.

vortex_detection.py

#   vortex_detection.py
#   (c) Daniel Bershatsky, 2016
#   See LICENSE for details.

from __future__ import absolute_import, division
import skimage.measure as measure
import numpy as np


def fit_segment(contour, index):
    A = np.array([contour[index + 0, :], contour[index + 1, :]])
    b = np.ones(2)
    return np.linalg.solve(A, b)

def check_intersection(first, i, second, j):
    v1, v2 = first[i, :], first[i + 1, :]
    u1, u2 = second[j, :], second[j + 1, :]

    a = v2 - v1
    b = u1 - v1
    c = u2 - v1
    d = u2 - u1

    return np.cross(a, b) * np.cross(a, c) <= 0.0 and np.cross(d, -b) * np.cross(d, -c) <= 0.0

def find_intercections(first, second, tol):
    intersections = list()

    for i in xrange(first.shape[0] - 1):
        for j in xrange(second.shape[0] - 1):
            if not check_intersection(first, i, second, j):
                continue

            u = fit_segment(first, i)
            v = fit_segment(second, j)

            A = np.array([u, v])
            b = np.ones(2)

            if np.linalg.matrix_rank(A) != 2:
                continue

            intersection = np.linalg.solve(A, b)
            angle = np.dot(u, v) / (np.linalg.norm(u) * np.linalg.norm(v))

            if any(np.linalg.norm(x[0] - intersection) < tol for x in intersections):
                continue

            intersections.append((intersection, angle))

    return intersections

def find_vortices_by_density(u, tol):
    reals = measure.find_contours(u.real, 0.0)
    imags = measure.find_contours(u.imag, 0.0)

    vortices = list()

    for i, real in enumerate(reals):
        for j, imag in enumerate(imags):
            vortices.extend(find_intercections(real, imag, tol))

    return zip(*vortices)

def find_vortices(u, tol=1.0e-6, by='density'):
    return find_vortices_by_density(u, tol)