"""Visualize that testing if a point p is to the left or right of a line
defined by two points q and r is nummerical unstable when using floats.

  q = (12, 12), r = (24, 24), p in [1/2, 1/2 + eps[^2 where eps = 2^-46

Attempt to reproduce slide 9 from
https://people.mpi-inf.mpg.de/~mehlhorn/ftp/SoCG09.pdf

Program plots all 720 possible permutations of adding the terms of the
determinant.
"""

import matplotlib.pyplot as plt
from matplotlib.backends.backend_pdf import PdfPages
from itertools import permutations

PDF_FILE = None
# PDF_FILE = "sign-plot.pdf"  # uncomment to generate PDF file

resolution = 5

N = 2**resolution
delta = 1 / 2**(46+resolution)

q = (12, 12)
r = (24, 24)

######################################################################
# Compute det values for all possible permutations of determinant terms

P = [N*N*[None] for _ in range(720)]  # P[permutation idx][point idx] = (x,y,z)

for i in range(N):
    for j in range(N):
        p = (0.5 + i*delta, 0.5 + j*delta)
        v = q[0]*r[1], r[0]*p[1], p[0]*q[1], -r[0]*q[1], -p[0]*r[1], -q[0]*p[1]

        for idx, v_permuted in enumerate(permutations(v)):
            P[idx][i*N+j] = (i, j, sum(v_permuted))

TERMS = ['q[0]*r[1]', 'r[0]*p[1]', 'p[0]*q[1]',
         '-r[0]*q[1]', '-p[0]*r[1]', '-q[0]*p[1]']

DETERMINANT_EQUATION = [' + '.join(v) for v in permutations(TERMS)]

######################################################################
# Plot data

cols = 12
rows = 6
subplts = cols * rows

if PDF_FILE:
    pdf = PdfPages(PDF_FILE)

for K in range(0, len(P), subplts):
    fig = plt.figure()
    plt.subplots_adjust(left=0.01, right=.99, bottom=0.01, top=0.99,
                        hspace=0.1, wspace=0.1)
    for k in range(K, min(len(P), K + subplts)):
        print(k, ":", DETERMINANT_EQUATION[k])

        ax = plt.subplot(rows, cols, 1 + k % subplts,
                         aspect='equal', facecolor='lightgrey')
        ax.axis([-0.5, N - 0.5, -0.5, N - 0.5]) 
        ax.set_title(str(k + 1), position=(.25, .65))
        ax.tick_params(bottom='off', labelbottom='off',
                       left='off', labelleft='off')

        positive = [(x, y) for x, y, z in P[k] if z > 0]
        negative = [(x, y) for x, y, z in P[k] if z < 0]
        zero = [(x, y) for x, y, z in P[k] if z == 0]

        for points, color in [(positive, "b"), (negative, "r"), (zero, "y")]:
            X = [x for x, y in points]
            Y = [y for x, y in points]
            plt.plot(X, Y, color + ".", ms=3)

        plt.plot([-1, N], [-1, N], "k-")

    if PDF_FILE:
        pdf.savefig()

if PDF_FILE:
    pdf.close()
else:
    plt.show()
