--- lua-tikz3dtools-vector.lua
--- Vector class for the lua-tikz3dtools package.

local Vector = {}
Vector.__index = Vector

--- Create a new 1xn Vector object (validated).
--- @param vec table<number> The 1xn vector data
--- @return Vector A 1xn vector
function Vector:new(vec)
    local A = type(vec) == "table"
    local n = #vec
    A = A and n > 0
    if A == true then
        for c = 1, n do
            A = A and type(vec[c]) == "number"
        end
    end
    assert(A, "A Vector object must be a 1xn table of numbers, where n>0.")
    return setmetatable(vec, Vector)
end

--- Unchecked constructor — skips validation. Use for internal calls only.
--- @param vec table<number>
--- @return Vector
function Vector:_new(vec)
    return setmetatable(vec, Vector)
end

--- Pretty print for Vector object
--- @return string
function Vector:__tostring()
    return ("Vector{%.12f" .. string.rep(",%.12f", #self - 1) .. "}"):format(table.unpack(self))
end

--- Deep copy of Vector object
--- @return Vector A deep copy of self
function Vector:deep_copy()
    local n = #self
    local copy = {}
    for i = 1, n do
        copy[i] = self[i]
    end
    return Vector:_new(copy)
end

--- Convert Vector to table
--- @return table<number> A table representation of self
function Vector:to_table()
    local n = #self
    local t = {}
    for i = 1, n do
        t[i] = self[i]
    end
    return t
end

--- Homogeneous addition of Vector objects
--- @param other Vector A vector of same size
--- @return Vector The sum of self and other
function Vector:hadd(other)
    local lo = #other
    local V = {}
    for i = 1, lo - 1 do
        V[i] = self[i] + other[i]
    end
    V[lo] = self[lo]
    return Vector:_new(V)
end

--- Addition of Vector objects
--- @param other Vector A vector of same size
--- @return Vector The sum of self and other
function Vector:add(other)
    local lo = #other
    local V = {}
    for i = 1, lo do
        V[i] = self[i] + other[i]
    end
    return Vector:_new(V)
end

--- Homogeneous subtraction of Vector objects
--- @param other Vector A vector of same size
--- @return Vector The difference of self and other
function Vector:hsub(other)
    local lo = #other
    local V = {}
    for i = 1, lo - 1 do
        V[i] = self[i] - other[i]
    end
    V[lo] = self[lo]
    return Vector:_new(V)
end

--- Subtraction of Vector objects
--- @param other Vector A vector of same size
--- @return Vector The difference of self and other
function Vector:sub(other)
    local lo = #other
    local V = {}
    for i = 1, lo do
        V[i] = self[i] - other[i]
    end
    return Vector:_new(V)
end

--- Homogeneous scaling of Vector objects
--- @param scalar number The scalar
--- @return Vector The scale of self by scalar
function Vector:hscale(scalar)
    local ls = #self
    local V = {}
    for i = 1, ls - 1 do
        V[i] = self[i] * scalar
    end
    V[ls] = self[ls]
    return Vector:_new(V)
end

--- Scaling of Vector objects
--- @param scalar number The scalar
--- @return Vector The scale of self by scalar
function Vector:scale(scalar)
    local ls = #self
    local V = {}
    for i = 1, ls do
        V[i] = self[i] * scalar
    end
    return Vector:_new(V)
end

--- Homogeneous inner product of Vector objects
--- @param other Vector The RHS
--- @return number The homogeneous inner product of self with other
function Vector:hinner(other)
    local ls = #self
    local result = 0
    for i = 1, ls - 1 do
        result = result + self[i] * other[i]
    end
    return result
end

--- Homogeneous hypercross product of Vector objects
--- @param ... Vector Additional vectors
--- @return Vector The homogeneous hypercross product of self with the additional vectors
function Vector:hhypercross(...)
    -- Forward-declared; Matrix is set via Vector._set_Matrix()
    local Matrix = Vector._Matrix
    local args = { ... }
    local n = #self - 1

    local M = {}
    for i = 1, n do
        local e = {}
        for j = 1, n do
            e[j] = (i == j) and 1 or 0
        end
        e[n + 1] = 0

        M[i] = {}
        M[i][1] = Vector:_new(e)
        M[i][2] = self[i]

        for j = 1, n - 2 do
            M[i][j + 2] = args[j][i]
        end
    end
    local result = Matrix:_new(M):det("column", 1)
    result[n + 1] = self[n + 1]
    return Vector:_new(result)
end

--- Homogeneous normalization of Vector objects
--- @return Vector The homogeneous normalization of self
function Vector:hnormalize()
    local ls = #self
    local len = 0
    for i = 1, ls - 1 do
        len = len + self[i] * self[i]
    end
    len = math.sqrt(len)
    local V = {}
    for i = 1, ls - 1 do
        V[i] = self[i] / len
    end
    V[ls] = self[ls]
    return Vector:_new(V)
end

--- Homogeneous distance between Vector objects
--- @param other Vector A vector of same size
--- @return number The homogeneous distance between self and other
function Vector:hdistance(other)
    local ls = #self
    local result = 0
    for i = 1, ls - 1 do
        result = result + (self[i] - other[i])^2
    end
    return math.sqrt(result)
end

--- Homogeneous projection of Vector objects
--- @param other Vector A vector of same size
--- @return Vector The homogeneous projection of self onto other
function Vector:hproject_onto(other)
    local denom = other:hinner(other)
    local scalar = self:hinner(other) / denom
    return other:hscale(scalar)
end

--- Homogeneous orthogonal projection of Vector objects onto a plane
--- @param plane_basis table A 3-row basis (origin, u-dir, v-dir) in homogeneous coordinates
--- @return Vector The homogeneous orthogonal projection of self onto the plane
function Vector:horthogonal_projection_onto_plane(plane_basis)
    local Matrix = Vector._Matrix
    local TO = Vector:_new(plane_basis[1])
    local TU = Vector:_new(plane_basis[2])
    local TV = Vector:_new(plane_basis[3])
    local TW = TU:hhypercross(TV)
    local TOS = self:hsub(TO)
    if TW:hinner(TOS) < 0 then
        TW = TW:hscale(-1)
    end
    local proj = TOS:hproject_onto(TW)
    return self:hsub(proj)
end

--- Homogeneous norm of Vector objects
--- @return number The homogeneous norm of self
function Vector:hnorm()
    local sum = 0
    for i = 1, #self - 1 do
        sum = sum + self[i] * self[i]
    end
    return math.sqrt(sum)
end

--- Homogeneous reciprocation of Vector objects
--- @return Vector The homogeneous reciprocation of self
function Vector:reciprocate_by_homogeneous()
    local V = {}
    local w = self[#self]
    for i = 1, #self - 1 do
        V[i] = self[i] / w
    end
    V[#self] = 1
    return Vector:_new(V)
end

--- Find an arbitrary vector orthogonal to self
--- @return Vector An arbitrary vector orthogonal to self
function Vector:orthogonal_vector()
    local v
    if (
        math.abs(self[2]) > 1e-3
        and math.abs(self[1]) < 1e-3
        and math.abs(self[3]) < 1e-3
    ) then
        v = self:hhypercross(Vector:_new{0,1,0,1})
    else
        v = self:hhypercross(Vector:_new{1,0,0,1})
    end
    return v
end

--- Multiply Vector by Matrix (applies projective divide)
--- @param matrix table The matrix
--- @return Vector The product of self and matrix, divided by homogeneous w
function Vector:multiply(matrix)
    local Matrix = Vector._Matrix
    local M = Matrix:_new{{table.unpack(self)}}:multiply(matrix)
    return Vector:_new(M[1]):reciprocate_by_homogeneous()
end

--- Create a homogeneous sphere point from spherical coordinates
--- @param v Vector the theta, phi, and radius
--- @return Vector
function Vector.sphere(v)
    local theta, phi, radius = v[1], v[2], v[3]
    local s = math.sin(phi)
    local x = radius * s * math.cos(theta)
    local y = radius * s * math.sin(theta)
    local z = radius * math.cos(phi)
    return Vector:_new{x, y, z, 1}
end

--- The stereographic projection
--- @param v Vector the 3D vector
--- @return Vector the stereographic projection
function Vector.stereographic_projection(v)
    local x, y, z = v[1], v[2], v[3]
    return Vector:_new{x/(1-z), y/(1-z), 0, 1}
end

--- The inverse stereographic projection
--- @param x number the x point
--- @param y number the y point
--- @return Vector the inverse stereographic projection
function Vector.inverse_stereographic_projection(v)
    local x, y = v[1], v[2]
    return Vector:_new{2*x/(1+x^2+y^2), 2*y/(1+x^2+y^2), (-1+x^2+y^2)/(1+x^2+y^2), 1}
end

--- Set the Matrix class reference (breaks circular dependency).
--- @param M table The Matrix class
function Vector._set_Matrix(M)
    Vector._Matrix = M
end

return Vector