#!/usr/bin/env python3
#
# line_equation - some experiments about representing a line
#
# Copyright (C) 2016  Antonio Ospite <ao2@ao2.it>
#
# 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 <http://www.gnu.org/licenses/>.

"""Some experiments about representing a line passing for two points"""

import math
import matplotlib.pyplot as plt
import numpy as np

EPSILON = 10e-6
INFINITY = 10e6

# pylint: disable=invalid-name


def distance(p1, p2):
    """Return the distance between the two points p1 and p2"""
    dx = p2[0] - p1[0]
    dy = p2[1] - p1[1]

    return math.sqrt(dx * dx + dy * dy)


def line_eq_slope_intercept(p1, p2):
    """Return the coefficients of the line equation in the slope/intercept
    form:

        y = mx + q
    """
    dx = p2[0] - p1[0]
    dy = p2[1] - p1[1]

    # For vertical lines the slope would be undefined in theory.
    # For practical purposes a large value can work.
    if dx == 0:
        slope = INFINITY * math.copysign(1, dy)
    else:
        slope = dy / dx

    intercept = -slope * p1[0] + p1[1]

    return slope, intercept


def line_angle(p1, p2):
    """Return the angle between the line passing for p1 and p2, and the positive
    x axis.

    The return value is in radians.
    """
    dx = p2[0] - p1[0]
    dy = p2[1] - p1[1]

    angle = math.atan2(dy, dx)
    if angle < 0:
        angle += 2 * math.pi

    return angle


def slope_intercept_to_a_b_c(slope, intercept):
    """Convert the coefficients from the slope/intercept form:

        y = mx + q

    to the coefficients in the linear equation form:

        ax + by + c = 0
    """
    a = -slope
    b = 1.
    c = -intercept

    return a, b, c


def line_eq_two_points(p1, p2):
    """Return the coefficients for the line equation, in the linear equation
    form:

        ax + by + c = 0
    """
    a = p2[1] - p1[1]
    b = - (p2[0] - p1[0])
    c = - (p1[0] * p2[1] - p2[0] * p1[1])

    return a, b, c


def line_calc_normal(p1, p2):
    """Return the unit normal to the line passing for two points."""
    dx = p2[0] - p1[0]
    dy = p2[1] - p1[1]

    n_length = math.sqrt(dx*dx + dy*dy)
    if n_length < EPSILON:
        nx = 0
        ny = 0
    else:
        nx = dy / n_length
        ny = -dx / n_length

    return nx, ny


def list_rotate(l, n):
    return l[n:] + l[:n]


def main():
    points = [
        [132.143, 941.429],
        [480.714, 905.714],
        [480.714, 735.714],
        [180.0, 783.571]
    ]

    # The points are assumed to be in clockwise order
    for p1, p2 in zip(points, list_rotate(points, 1)):
        print("p1", p1)
        print("p2", p2)

        line_x_extremes = [p1[0], p2[0]]
        line_y_extremes = [p1[1], p2[1]]

        plt.plot(line_x_extremes, line_y_extremes, 'k')

        x_samples = np.linspace(min(line_x_extremes), max(line_x_extremes), 15)

        # slope/intercept form
        m, q = line_eq_slope_intercept(p1, p2)

        if abs(m) == INFINITY:
            print("m is undefined: vertical line!")
            y_samples = \
                np.linspace(min(line_y_extremes), max(line_y_extremes), 15)
        else:
            y_samples = x_samples * m + q

        plt.plot(x_samples, y_samples, 'wo', markersize=15)

        # linear equation form
        a, b, c = line_eq_two_points(p1, p2)

        if b == 0:
            print("b is 0: vertical line!")
            y_samples = \
                np.linspace(min(line_y_extremes), max(line_y_extremes), 15)
        else:
            y_samples = (-x_samples * a - c) / b

        plt.plot(x_samples, y_samples, 'ws', markersize=10)

        # parametric forms
        t = np.linspace(0, 1, len(x_samples))

        # point and direction form
        dx = p2[0] - p1[0]
        dy = p2[1] - p1[1]

        x_samples = p1[0] + dx * t
        y_samples = p1[1] + dy * t
        plt.plot(x_samples, y_samples, 'm+', markersize=10)

        # or using a point and the normal, but in this case the normal has to
        # be rotated by 90 degrees counterclockwise
        nx, ny = line_calc_normal(p1, p2)
        rnx = -ny
        rny = nx

        length = distance(p1, p2)
        x_samples = p1[0] + rnx * t * length
        y_samples = p1[1] + rny * t * length
        plt.plot(x_samples, y_samples, 'bx', markersize=10)

        # draw the normal
        mid_x = (p1[0] + p2[0]) / 2
        mid_y = (p1[1] + p2[1]) / 2
        n_length = 20
        plt.arrow(mid_x, mid_y, nx * n_length, ny * n_length,
                  head_width=10, head_length=5, fc='k', ec='k')

        # check the other functions
        print("angle: %g" % math.degrees(line_angle(p1, p2)))
        a1, b1, c1 = slope_intercept_to_a_b_c(m, q)
        print(a1, b1, c1)
        if b != 0:
            print(a / b, b / b, c / b)
        else:
            print(a, b, c)
        print()

    # draw the vertices on top of everything else
    xvalues = [p[0] for p in points]
    yvalues = [p[1] for p in points]
    plt.plot(xvalues, yvalues, 'ro')
    plt.title('Line passing for two points')
    plt.xlabel('x')
    plt.ylabel('y')

    plt.plot([], 'k', label="line")
    plt.plot([], 'wo', label="y = mx + q")
    plt.plot([], 'ws', label="ax + by + c = 0")
    plt.plot([], 'm+', label="P = P1 + D * t")
    plt.plot([], 'bx', label="P = P1 + N * t * length")
    plt.legend(loc='lower left', numpoints=1, fancybox=True, prop={'size': 10})
    plt.show()


if __name__ == "__main__":
    main()