using Isoband
using Test

# let's use this quadruple bar thing as the replacement for setequal() from the R tests for convenience
≣(s1::Set, s2::Set) = s1 == s2
≣(s1, s2) = Set(s1) == Set(s2)

@testset "elementary polygons get merged" begin
    # two connected polygons get merged
    z = [0 0 1;
         1 1 1]

    out = isobands(1:3, 2:-1:1, z, 0.5, 1.5)

    @test 1000 * out.x + out.y ≣ 1000 * [3.0, 2.0, 1.0, 1.0, 2.0, 2.5, 3.0] + [1.0, 1.0, 1.0, 1.5, 1.5, 2.0, 2.0]
    @test out.id == ones(7)

    # two unconnected polygons don't get merged
    z = [1 2 1;
         1 2 2]
    out = isobands(1:3, 2:-1:1, z, 0.5, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [3.0, 2.5, 3.0, 1.0, 1.5, 1.5, 1.0] + [1.5, 2.0, 2.0, 2.0, 2.0, 1.0, 1.0]
    @test out.id ≣ [1, 2]
    @test length(out.id) == 7

    # two separate bands get merged in second row
    z = [1 2 1;
         1 2 1;
         0 0 0]

    out = isobands(1:3, 3:-1:1, z, 0.5, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [3.0, 2.0, 1.0, 1.0, 1.0, 1.5, 1.5, 2.0, 2.5, 2.5, 3.0, 3.0] + [1.50, 1.25, 1.50, 2.00, 3.00, 3.00, 2.00, 1.75, 2.00, 3.00, 3.00, 2.00]
    @test out.id == ones(12)

    # circle gets closed
    z = [1 1 1;
         1 2 1;
         1 1 1]

    out = isobands(1:3, 3:-1:1, z, 0.5, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [3.0, 2.0, 1.0, 1.0, 1.0, 2.0, 3.0, 3.0, 2.0, 2.5, 2.0, 1.5] + [1.0, 1.0, 1.0, 2.0, 3.0, 3.0, 3.0, 2.0, 1.5, 2.0, 2.5, 2.0]
    @test out.id ≣ [1, 2]
    @test length(out.id) == 12
end

@testset "NaNs are handled correctly" begin
    z = [NaN 1 1;
         1 1 1;
         1 1 1]

    out = isobands(1:3, 3:-1:1, z, 0.5, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [3.0, 2, 1, 1, 2, 2, 3, 3] + [1.0, 1, 1, 2, 2, 3, 3, 2]
    @test out.id == ones(8)
  
    z = [NaN 1 1;
         1 1 1;
         1 1 NaN]

    out = isobands(1:3, 3:-1:1, z, 0.5, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [1.0, 1, 2, 2, 2, 2, 3, 3] + [1, 2, 2, 1, 2, 3, 3, 2]
    @test out.id ≣ [1, 2]
    @test length(out.id) == 8
end

@testset "All elementary shapes are calculated correctly" begin
  # a matrix that requires all elementary shapes for isobanding
  z =  [0 0 0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2;
        0 0 0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2;
        0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2 0 0;
        1 0 1 1 1 2 2 0 2 1 2 2 0 0 0 1 0 2;
        2 0 2 1 2 2 0 0 0 1 0 2 1 0 1 1 1 2;
        0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2 0 0;
        2 0 2 1 2 2 0 0 0 1 0 2 1 0 1 1 1 2;
        1 0 1 1 1 2 2 0 2 1 2 2 0 0 0 1 0 2;
        0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2 0 0;
        0 0 0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2;
        0 0 0 1 0 2 1 0 1 1 1 2 2 0 2 1 2 2]

    out = isobands(1:18, 11:-1:1, z, 0.5, 1.5)

    @test 1000 * out.x + out.y ≣ 1000 * [
        16.00, 15.50, 15.50, 16.00, 16.25, 16.00, 15.00, 14.50, 14.75, 14.75, 14.25, 14.25,
        14.00, 13.75, 13.75, 13.25, 13.25, 13.50, 13.00, 12.75, 12.00, 11.50, 11.50, 11.00,
        10.00,  9.00,  8.50,  8.50,  8.00,  7.50,  7.50,  7.00,  6.50,  6.50,  6.00,  5.75,
        5.75,  5.25,  5.25,  5.00,  4.50,  4.50,  4.00,  3.50,  3.50,  3.25,  3.00,  2.50,
        2.00,  1.50,  1.00,  1.00,  1.00,  1.25,  1.00,  1.00,  1.50,  1.00,  1.00,  1.25,
        1.00,  1.00,  1.00,  1.50,  2.00,  2.50,  3.00,  3.25,  3.50,  3.50,  4.00,  4.50,
        4.50,  5.00,  5.25,  5.25,  5.75,  5.75,  6.00,  6.50,  6.50,  7.00,  7.50,  7.50,
        8.00,  8.50,  8.50,  9.00, 10.00, 11.00, 11.50, 11.50, 12.00, 12.75, 13.00, 13.50,
        13.25, 13.25, 13.75, 13.75, 14.00, 14.25, 14.25, 14.75, 14.75, 14.50, 15.00, 16.00,
        16.25, 16.00, 15.50, 15.50, 16.00, 16.50, 16.50, 17.00, 18.00, 18.00, 17.00, 16.75,
        16.50, 17.00, 17.25, 18.00, 18.00, 17.75, 17.50, 18.00, 18.00, 17.00, 16.75, 17.00,
        18.00, 18.00, 17.50, 17.75, 18.00, 18.00, 17.25, 17.00, 16.50, 16.75, 17.00, 18.00,
        18.00, 17.00, 16.50, 16.50, 12.00, 12.25, 12.00, 11.75, 11.00, 11.25, 12.00, 12.25,
        12.50, 12.00, 11.75, 11.00, 10.50, 10.00,  9.50,  9.00,  8.75,  9.00,  9.50, 10.00,
        6.00,  6.50,  6.00,  5.50,  4.00,  4.50,  4.00,  3.75, 13.00, 14.00, 15.00, 15.50,
        15.00, 14.50, 14.00, 13.50, 13.00, 12.75,  8.00,  8.25,  9.00,  9.50,  9.00,  8.00,
        7.00,  6.75,  7.00,  7.75,  7.00,  7.25,  7.00,  6.25,  6.00,  5.00,  4.50,  5.00,
        6.00,  6.25,  7.00,  7.25,  7.00,  6.00,  5.50,  5.00,  4.50,  4.00,  3.75,  4.00,
        4.50,  5.00,  5.50,  6.00,  2.00,  2.50,  2.25,  2.00,  1.75,  1.50, 11.00, 11.25,
        11.00, 10.50, 13.00, 13.50, 13.00, 12.75, 12.00, 12.25, 12.00, 11.75, 11.00, 11.25,
        11.00, 10.00,  9.50, 10.00,  9.00,  8.25,  8.00,  7.75,  7.00,  6.75,  7.00,  8.00,
        9.00,  9.50, 16.25, 16.00, 15.00, 14.50, 15.00, 16.00, 15.00, 14.00, 13.00, 12.75,
        13.00, 13.50, 14.00, 14.50, 15.00, 15.50, 12.00, 12.50, 12.25, 12.00, 11.25, 11.00,
        10.00,  9.50,  9.00,  8.75,  9.00,  9.50, 10.00, 10.50, 11.00, 11.75, 11.00, 11.25,
        11.00, 10.50,  6.00,  5.50,  6.00,  6.50,  3.00,  2.50,  3.00,  3.25,  3.00,  3.50,
        3.00,  2.75,  2.00,  2.25,  2.50,  2.00,  1.50,  1.75,  3.00,  2.75,  3.00,  3.50,
        4.00,  3.75,  4.00,  4.50, 12.00, 12.25, 12.00, 11.75,  6.00,  6.50,  6.00,  5.50] +
        
        [
        1.00,  1.00,  2.00,  2.50,  3.00,  3.50,  3.25,  3.00,  2.00,  1.00,  1.00,  2.00,
        2.50,  2.00,  1.00,  1.00,  2.00,  3.00,  3.25,  3.00,  2.25,  2.00,  1.00,  1.00,
        1.00,  1.00,  1.00,  2.00,  2.50,  2.00,  1.00,  1.00,  1.00,  2.00,  2.25,  2.00,
        1.00,  1.00,  2.00,  2.50,  2.00,  1.00,  1.00,  1.00,  2.00,  3.00,  3.50,  3.00,
        2.50,  3.00,  3.50,  4.00,  4.50,  5.00,  5.25,  5.75,  6.00,  6.25,  6.75,  7.00,
        7.50,  8.00,  8.50,  9.00,  9.50,  9.00,  8.50,  9.00, 10.00, 11.00, 11.00, 11.00,
       10.00,  9.50, 10.00, 11.00, 11.00, 10.00,  9.75, 10.00, 11.00, 11.00, 11.00, 10.00,
        9.50, 10.00, 11.00, 11.00, 11.00, 11.00, 11.00, 10.00,  9.75,  9.00,  8.75,  9.00,
       10.00, 11.00, 11.00, 10.00,  9.50, 10.00, 11.00, 11.00, 10.00,  9.00,  8.75,  8.50,
        9.00,  9.50, 10.00, 11.00, 11.00, 11.00, 10.00,  9.75,  9.75,  9.25,  9.25,  9.00,
        8.00,  7.50,  8.00,  8.75,  8.25,  8.00,  7.00,  6.75,  6.25,  6.50,  6.00,  5.50,
        5.75,  5.25,  5.00,  4.00,  3.75,  3.25,  4.00,  4.50,  4.00,  3.00,  2.75,  2.75,
        2.25,  2.25,  2.00,  1.00,  2.75,  3.00,  3.25,  3.00,  2.50,  3.00,  3.75,  4.00,
        5.00,  5.25,  5.00,  4.25,  4.00,  3.50,  4.00,  4.25,  4.00,  3.50,  3.00,  2.50,
        2.75,  3.00,  3.25,  3.00,  2.50,  3.00,  3.50,  3.00,  3.75,  3.50,  3.75,  4.00,
        4.50,  5.00,  5.50,  5.00,  4.50,  4.00,  3.50,  4.00,  4.75,  5.00,  5.50,  5.50,
        5.50,  5.00,  4.75,  4.00,  3.50,  4.00,  4.25,  5.00,  5.25,  5.50,  6.00,  6.50,
        6.75,  7.00,  7.75,  8.00,  8.50,  8.25,  8.00,  7.50,  7.00,  6.50,  6.00,  5.50,
        5.00,  4.50,  4.00,  3.75,  3.50,  4.00,  5.00,  5.50,  5.00,  4.00,  4.75,  5.00,
        5.25,  5.00,  5.50,  6.00,  6.50,  6.00,  5.75,  6.00,  6.25,  6.00,  5.75,  6.00,
        6.25,  6.50,  6.00,  5.50,  7.25,  8.00,  8.50,  8.00,  7.25,  7.00,  6.50,  6.50,
        6.50,  7.00,  6.00,  6.50,  6.50,  6.00,  5.50,  5.50,  8.25,  8.50,  8.25,  8.00,
        7.50,  7.00,  6.50,  7.00,  7.50,  8.00,  6.75,  7.00,  8.00,  8.25,  9.00,  9.50,
        9.50,  9.00,  8.50,  8.00,  7.75,  8.00,  8.50,  8.00,  7.75,  7.00,  6.75,  7.00,
        7.25,  7.00,  6.25,  6.00,  5.75,  6.00,  6.25,  6.00,  5.75,  6.00,  6.75,  7.00,
        7.50,  7.00,  6.50,  7.00,  8.00,  8.50,  8.00,  7.00,  5.25,  5.00,  4.50,  5.00,
        9.50,  9.00,  8.50,  9.00,  8.75,  9.00,  9.25,  9.00,  8.75,  9.00,  9.25,  9.00
     ]
  
  @test out.id ≣ 1:26
  @test length(out.id) == 324
end


@testset "Six-sided saddles" begin
    # a matrix that contains all six-sided saddles
    z = [0 1 1 2;
        1 0 2 1;
        0 1 1 2]
        
    # midpoint outside the band
    out = isobands(1:4, 3:-1:1, z, 0.6, 1.4)
    @test 1000 * out.x + out.y ≣ 1000 * [3.0, 2.0, 1.6, 2.0, 2.3, 2.0,
        1.6, 2.0, 3.0, 3.4, 3.0, 2.7,
        3.0, 3.4, 1.0, 1.0, 1.0, 1.4,
        4.0, 4.0, 4.0, 3.6] + 
     [1.0, 1.0, 1.0, 1.4, 2.0, 2.6,
        3.0, 3.0, 3.0, 3.0, 2.6, 2.0,
        1.4, 1.0, 1.6, 2.0, 2.4, 2.0,
        2.4, 2.0, 1.6, 2.0]
    @test out.id ≣ 1:3
    @test length(out.id) == 22
  
    # midpoint inside the band
    out = isobands(1:4, 3:-1:1, z, 0.4, 1.6)
    @test 1000 * out.x + out.y ≣ 1000 * [3.0, 2.0, 1.4, 1.0, 1.0, 1.0,
        1.4, 2.0, 3.0, 3.6, 4.0, 4.0,
        4.0, 3.6, 3.0, 3.4, 3.0, 2.8,
        2.0, 2.2, 2.0, 1.6] +
    [1.0, 1.0, 1.0, 1.4, 2.0, 2.6,
        3.0, 3.0, 3.0, 3.0, 2.6, 2.0,
        1.4, 1.0, 1.6, 2.0, 2.4, 2.0,
        1.6, 2.0, 2.4, 2.0] 
    @test out.id ≣ 1:3
    @test length(out.id) == 22
end


@testset "Seven-sided saddles" begin
    # a matrix that contains all seven-sided saddles
    z = [0 1 0 1 2 1;
        2 0 2 2 0 2;
        0 1 0 1 2 1]

    # midpoint inside the band
    out = isobands(1:6, 3:-1:1, z, 0.5, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [6.00, 6.00, 5.50, 5.00, 4.50, 4.00,
        3.50, 3.00, 2.50, 2.00, 1.50, 1.00,
        1.00, 1.25, 1.00, 1.00, 1.50, 2.00,
        2.50, 3.00, 3.50, 4.00, 4.50, 5.00,
        5.50, 6.00, 6.00, 5.75, 4.00, 4.25,
        4.00, 3.00, 2.75, 3.00, 2.00, 2.25,
        2.00, 1.75, 5.00, 5.25, 5.00, 4.75] +
    [1.50, 1.00, 1.00, 1.25, 1.00, 1.00,
        1.00, 1.25, 1.00, 1.00, 1.00, 1.25,
        1.75, 2.00, 2.25, 2.75, 3.00, 3.00,
        3.00, 2.75, 3.00, 3.00, 3.00, 2.75,
        3.00, 3.00, 2.50, 2.00, 1.50, 2.00,
        2.50, 2.25, 2.00, 1.75, 1.50, 2.00,
        2.50, 2.00, 1.75, 2.00, 2.25, 2.00]
    
    @test out.id ≣ 1:4
    @test length(out.id) == 42
  
    # midpoint outside the band
    out = isobands(1:6, 3:-1:1, z, 0.8, 1.2)
    @test 1000 * out.x + out.y ≣ 1000 * [6.0, 6.0, 5.8, 5.0, 4.4, 5.0,
        5.6, 4.2, 4.0, 3.8, 3.0, 2.4,
        3.0, 3.8, 4.0, 4.2, 4.0, 3.0,
        2.6, 3.0, 4.0, 2.2, 2.0, 1.8,
        2.0, 1.0, 1.4, 1.0, 1.0, 1.6,
        1.0, 5.0, 5.4, 5.0, 4.6, 6.0,
        5.8, 6.0, 2.2, 2.0, 1.8, 2.0] +
    [1.2, 1.0, 1.0, 1.4, 2.0, 2.6,
        2.0, 1.0, 1.0, 1.0, 1.4, 2.0,
        2.6, 3.0, 3.0, 3.0, 2.8, 2.4,
        2.0, 1.6, 1.2, 1.0, 1.0, 1.0,
        1.2, 1.6, 2.0, 2.4, 2.6, 2.0,
        1.4, 1.6, 2.0, 2.4, 2.0, 2.8,
        3.0, 3.0, 3.0, 2.8, 3.0, 3.0]
    @test out.id ≣ 1:8
    @test length(out.id) == 42
end

@testset "Eight-sided saddles" begin
    # a matrix that contains all eight-sided saddles
    z = [0 2 0;
         2 0 2]

    # midpoint above the band
    out = isobands(1:3, 2:-1:1, z, 0.5, 0.8)
    @test 1000 * out.x + out.y ≣ 1000 * [3.00, 3.00, 2.60, 2.75, 2.00,
        2.40, 2.25, 2.00, 1.75, 1.60,
        1.00, 1.00, 1.25, 1.40] +
    [1.75, 1.60, 2.00, 2.00, 1.40,
        1.00, 1.00, 1.25, 1.00, 1.00,
        1.60, 1.75, 2.00, 2.00]
    @test out.id ≣ 1:3
    @test length(out.id) == 14
  
    # midpoint inside the band
    out = isobands(1:3, 2:-1:1, z, 0.5, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [2.75, 2.25, 2.00, 1.75, 1.25,
        1.00, 1.00, 1.25, 1.75, 2.00,
        2.25, 2.75, 3.00, 3.00] +
    [1.00, 1.00, 1.25, 1.00, 1.00,
        1.25, 1.75, 2.00, 2.00, 1.75,
        2.00, 2.00, 1.75, 1.25]
    @test out.id == ones(14)
  
    # midpoint below the band
    out = isobands(1:3, 2:-1:1, z, 1.2, 1.5)
    @test 1000 * out.x + out.y ≣ 1000 * [3.00, 3.00, 2.75, 2.60, 2.40,
        2.00, 1.60, 1.75, 2.00, 2.25,
        1.40, 1.25, 1.00, 1.00] +
    [1.40, 1.25, 1.00, 1.00, 2.00,
        1.60, 2.00, 2.00, 1.75, 2.00,
        1.00, 1.00, 1.25, 1.40]
    @test out.id ≣ 1:3
    @test length(out.id) == 14
end

@testset "Inconsistent numbers of isoband levels cause an error" begin
  m = [0 0 1 1;
        0 1 1 1;
        1 1 0 0;
        0 0 1 0]

  # error, multiple values are not recycled
  @test_throws ErrorException isobands(1:4, 1:4, m, [0.5, 1.5, 2.5], [0.5, 1.5])
  @test_throws ErrorException isobands(1:4, 1:4, m, [0.5, 1.5], [0.5, 1.5, 2.5])
end



#######################################################
# isolines

@testset "line segments get merged" begin
    # two connected line segments get merged
    z = [0 0 1;
        1 1 1]
    out = isolines(1:3, 2:-1:1, z, 0.5)
    @test 1000 * out.x + out.y ≣ 1000 * [1, 2, 2.5] + [1.5, 1.5, 2.0]
    @test out.id == ones(3)

    # two unconnected line segments don't get merged
    z = [0 1 0;
        0 1 1]
    out = isolines(1:3, 2:-1:1, z, 0.5)
    @test 1000 * out.x + out.y ≣ 1000 * [2.5, 3.0, 1.5, 1.5] + [2.0, 1.5, 2.0, 1.0]
    @test out.id ≣ 1:2
    @test length(out.id) == 4

    # two separate lines get merged in second row
    z = [0 1 0;
        0 1 0;
        0 0 0]
    out = isolines(1:3, 3:-1:1, z, 0.5)
    @test 1000 * out.x + out.y ≣ 1000 * [2.5, 2.5, 2.0, 1.5, 1.5] + [3.0, 2.0, 1.5, 2.0, 3.0]
    @test out.id == ones(5)

    # circle gets closed
    z = [0 0 0;
        0 1 0;
        0 0 0]
    out = isolines(1:3, 3:-1:1, z, 0.5)
    # circle is closed
    @test out.x[1] == out.x[5]
    @test out.y[1] == out.y[5]
    # coords are correct
    @test 1000 * out.x + out.y ≣ 1000 * [2.5, 2.0, 1.5, 2.0, 2.5] + [2.0, 2.5, 2.0, 1.5, 2.0]
    @test out.id == ones(5)
end

@testset "NAs are handled correctly" begin
    z = [NaN 0 0;
        0 1 1;
        0 1 1]
    out = isolines(1:3, 3:-1:1, z, 0.5)
    @test 1000 * out.x + out.y ≣ 1000 * [1.5, 1.5, 2.0, 3.0] + [2.0, 1.0, 2.5, 2.5]
    @test out.id ≣ 1:2
    @test length(out.id) == 4
end

@testset "All elementary segments are calculated correctly" begin
    # a matrix that requires all elementary segments for isolines
    z = [0 0 0 1 1 0 1 1;
        0 0 0 1 1 0 1 1;
        0 1 1 0 1 1 0 0;
        1 1 0 0 0 1 1 0;
        1 0 1 1 0 0 0 1]
    out = isolines(1:8, 5:-1:1, z, 0.5)

    @test 1000 * out.x + out.y ≣ 1000 * [7.5, 7.0, 6.0, 5.5, 5.0, 4.5, 4.0, 3.5,
            3.0, 2.5, 3.0, 4.0, 4.5, 2.5, 2.0, 1.5,
            8.0, 7.0, 6.5, 7.0, 7.5, 8.0, 6.5, 6.5,
            6.0, 5.5, 5.5, 3.5, 3.5, 3.0, 2.0, 1.5,
            1.0] +
            [1.0, 1.5, 1.5, 2.0, 2.5, 3.0, 3.5, 3.0,
            2.5, 2.0, 1.5, 1.5, 1.0, 1.0, 1.5, 1.0,
            3.5, 3.5, 3.0, 2.5, 2.0, 1.5, 5.0, 4.0,
            3.5, 4.0, 5.0, 5.0, 4.0, 3.5, 3.5, 3.0,
            2.5]

    @test out.id ≣ 1:5
    @test length(out.id) == 33
end


@testset "Saddles" begin
  # a matrix that contains all saddles (there are only two)
  z = [0 1 0;
       1 0 1]

  out = isolines(1:3, 2:-1:1, z, 0.5)

  @test 1000 * out.x + out.y ≣ 1000 * [2.5, 3.0, 2.5, 2.0, 1.5, 1.5, 1.0] + [2.0, 1.5, 1.0, 1.5, 1.0, 2.0, 1.5]
  @test out.id ≣ 1:3
  @test length(out.id) == 7

  out = isolines(1:3, 2:-1:1, z, 0.6)

  @test 1000 * out.x + out.y ≣ 1000 * [1.6, 2.0, 2.4, 3.0, 2.6, 1.0, 1.4] + [2.0, 1.6, 2.0, 1.4, 1.0, 1.4, 1.0]
  @test out.id ≣ 1:3
  @test length(out.id) == 7
end