Skip to content
Draft
40 changes: 27 additions & 13 deletions src/Contour.jl
Original file line number Diff line number Diff line change
Expand Up @@ -66,6 +66,13 @@ function contour(x, y, z, level::Number)
trace_contour(x, y, z, level, cells, cell_pop)
end

function contour(cells, x, y, z, level::Number)
# Todo: size checking on x,y,z
cell_pop = get_level_cells!(cells, z, level)
trace_contour(x, y, z, level, cells, cell_pop)
end


"""
`contours` returns a set of isolines.

Expand All @@ -77,7 +84,10 @@ contours(::Any...) = error("This method exists only for documentation purposes")
`contours(x,y,z,levels)` Trace the contour levels indicated by the `levels`
argument.
"""
contours(x, y, z, levels) = ContourCollection([contour(x, y, z, l) for l in levels])
function contours(x, y, z, levels)
cells = cell_matrix(z)
ContourCollection([contour(cells, x, y, z, l) for l in levels])
end

"""
`contours(x,y,z,Nlevels::Integer)` Trace `Nlevels` contour levels at heights
Expand Down Expand Up @@ -187,6 +197,7 @@ end

# Maps cell type to crossing types for non-ambiguous cells
const edge_LUT = (SW, SE, EW, NE, 0x0, NS, NW, NW, NS, 0x0, NE, EW, SE, SW)
const start_edge_LUT = (0x00, 0x00, 0x01, 0x00, 0x01, 0x02, 0x00, 0x00, 0x01, 0x02, 0x00, 0x04, 0x00, 0x00)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think these should probably be refactored to use some name constants instead of just magic numbers.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for looking at this. I agree. Even better, on the latest diff these aren't there at all.

I was playing around with collapsing the branching logic but it turns out that the codegen was the same since Julia inlined and LLVM unrolled the loops and merged everything.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yeah, I haven't actually spent any time with this codebase for several years, so someone else should probably review it too. That's why I left it as a flat comment, without either 👍 or 👎 - I want to leave the actual decision to someone else... :)


function _get_case(z, h)
case = z[1] > h ? 0x01 : 0x00
Expand All @@ -196,10 +207,19 @@ function _get_case(z, h)
case
end

function get_level_cells(z, h::Number)

function cell_matrix(z)
xi_max, yi_max = size(z)
cells = zeros(UInt8, (xi_max-1,yi_max-1))
end

function get_level_cells(z,h::Number)
cells = cell_matrix(z)
cell_pop = get_level_cells!(cells, z, h)
cells, cell_pop
end

function get_level_cells!(cells, z, h::Number)
xi_max, yi_max = size(z)

cell_pop = 0

Expand Down Expand Up @@ -235,12 +255,12 @@ function get_level_cells(z, h::Number)
end
end

return cells, cell_pop
return cell_pop
end

function findfirst_cell(m, from_x, from_y)
s = size(m)
for xi = from_x:s[1], yi = from_y:s[2]
@inbounds for xi = from_x:s[1], yi = from_y:s[2]
!iszero(m[xi,yi]) && return xi,yi
end
return 0,0
Expand Down Expand Up @@ -354,7 +374,7 @@ function trace_contour(x, y, z, h::Number, cells::Array, cell_pop)
nonempty_cells = cell_pop
(xi_0, yi_0) = (1,1)

while nonempty_cells > 0
@inbounds while nonempty_cells > 0
contour = Curve2(promote_type(map(eltype, (x, y, z))...))

# Pick initial box
Expand All @@ -365,13 +385,7 @@ function trace_contour(x, y, z, h::Number, cells::Array, cell_pop)

# Pick a starting edge
crossing = get_first_crossing(cell)
starting_edge = UInt8(0)
for edge in (N, S, E, W)
if !iszero(edge & crossing)
starting_edge = edge
break
end
end
starting_edge = start_edge_LUT[crossing]

# Add the contour entry location for cell (xi_0,yi_0)
add_vertex!(contour, interpolate(x, y, z, h, xi_0, yi_0, starting_edge), fwd)
Expand Down