diff --git a/LICENSE.md b/LICENSE.md index 76745a6..f4b6aef 100644 --- a/LICENSE.md +++ b/LICENSE.md @@ -21,14 +21,13 @@ The Devices.jl package is licensed under the MIT "Expat" License: > TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE > SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. -A modified version of `adaptive_grid` was originally from PlotUtils.jl, -an MIT "Expat" licensed Julia package by Tom Breloff and contributors (notably -Kristoffer Carlsson). +A modified version of `adaptive_grid` was originally from PlotUtils.jl, an MIT "Expat" +licensed Julia package by Tom Breloff and contributors (notably Kristoffer Carlsson). -Cadence Design Systems, Inc. holds the rights to the GDS-II format. The specification -has been described with permission in the SPIE Handbook of Microlithography, Micromachining -and Microfabrication, vol. 1 (accessible [here](http://www.cnf.cornell.edu/cnf_spie9.html) -as of April 17, 2017). +Cadence Design Systems, Inc. holds the rights to the GDS-II format. The specification has +been described with permission in the SPIE Handbook of Microlithography, Micromachining and +Microfabrication, vol. 1 (accessible [here](http://www.cnf.cornell.edu/cnf_spie9.html) as of +April 17, 2017). This package uses the D3.js library for interactive visualization, under the 3-clause BSD license: @@ -58,3 +57,6 @@ license: > ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT > (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS > SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + +The file `src/predicates.jl` contains a Julia adaptation of public domain C code written by +Jonathan Richard Shewchuk. Available at http://www.cs.cmu.edu/~quake/robust.html. diff --git a/NEWS.md b/NEWS.md index 6bf0693..ffeb292 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,14 @@ +- v0.5.0 + - Rewrite some line/ray/line-segment intersection code. + - `adjust!` has been renamed `reconcile!` + - For many array-like methods on paths, there is now a keyword argument that lets you + specify whether you want to reconcile the path immediately or defer. This is useful + when speed is a concern and you are chaining operations. + - Implement code to automatically handle path intersections: `Paths.intersect!` + - Path termination: `Paths.terminate!` + - Style translation: `Paths.pin` and `Paths.translate`. + - Node/segment/style splitting: `Paths.split`. + - Add some default colors: 21--29 is a red gradient, 31--39 is green, 41--49 is blue. - v0.4.0 - Julia 1.0 compatibility. - `Rectangles.isproper` has a slightly different definition. Rectangles are diff --git a/REQUIRE b/REQUIRE index f02b1ad..709df52 100644 --- a/REQUIRE +++ b/REQUIRE @@ -7,3 +7,5 @@ Clipper 0.5.0 Unitful 0.11.0 Cairo 0.5.4 WebIO 0.4 0.5- +IntervalSets 0.3.1 +DataStructures 0.14 diff --git a/src/Devices.jl b/src/Devices.jl index eb2c9c0..dbe4c29 100644 --- a/src/Devices.jl +++ b/src/Devices.jl @@ -1,5 +1,6 @@ module Devices using Random +using LinearAlgebra using ForwardDiff using FileIO using WebIO @@ -10,8 +11,9 @@ import Clipper import Clipper: cclipper import FileIO: save, load -import Base: length, show, eltype -import Unitful: Length, DimensionlessQuantity, ustrip, unit, inch +import Base: length, show, eltype, intersect! +import Unitful: Length, LengthUnits, DimensionlessQuantity, NoUnits, DimensionError +import Unitful: ustrip, unit, inch Unitful.@derived_dimension InverseLength inv(Unitful.𝐋) export GDSMeta @@ -21,6 +23,31 @@ const DEFAULT_LAYER = 0 const DEFAULT_DATATYPE = 0 const D3STR = read(joinpath(dirname(@__FILE__), "../deps/d3.min.js"), String) +# setup for robust 2d predicates +const splitter, epsilon = + let every_other = true, half = 0.5, epsilon = 1.0, splitter = 1.0, check = 1.0 + lastcheck = check + epsilon *= half + every_other && (splitter *= 2.0) + every_other = !every_other + check = 1.0 + epsilon + + while (check != 1.0) && (check != lastcheck) + lastcheck = check + epsilon *= half + every_other && (splitter *= 2.0) + every_other = !every_other + check = 1.0 + epsilon + end + splitter += 1.0 + splitter, epsilon + end + +const resulterrbound = (3.0 + 8.0 * epsilon) * epsilon +const ccwerrboundA = (3.0 + 16.0 * epsilon) * epsilon +const ccwerrboundB = (2.0 + 12.0 * epsilon) * epsilon +const ccwerrboundC = (9.0 + 64.0 * epsilon) * epsilon * epsilon + function __init__() # To ensure no crashes global _clip = Ref(Clipper.Clip()) @@ -51,6 +78,12 @@ Type alias for numeric types suitable for coordinate systems. """ const Coordinate = Union{Real, Length} +""" + CoordinateUnits = Union{typeof(NoUnits), LengthUnits} +Type alias for units suitable for coordinate systems. +""" +const CoordinateUnits = Union{typeof(NoUnits), LengthUnits} + """ PointTypes = Union{Real, DimensionlessQuantity, Length, InverseLength} Allowed type variables for `Point{T}` types. @@ -85,6 +118,8 @@ export Points, Point, Rotation, Translation, XReflection, YReflection, gety, ∘ +include("predicates.jl") + abstract type Meta end struct GDSMeta <: Meta layer::Int @@ -146,6 +181,10 @@ export Cells, Cell, CellArray, CellPolygon, CellReference, uniquename include("utils.jl") + +function bridge! end # needed in Paths, but not defined in Microwave +export bridge! + include("paths/paths.jl") import .Paths: Path, α0, @@ -174,6 +213,7 @@ import .Paths: Path, style0, style1, setstyle!, + terminate!, turn!, undecorated @@ -200,6 +240,7 @@ export Paths, Path, Segment, Style, style1, setstyle!, straight!, + terminate!, turn!, undecorated @@ -207,7 +248,6 @@ include("render/render.jl") include("microwave.jl") import .Microwave: - bridge!, checkerboard!, device_template!, flux_bias!, @@ -221,7 +261,6 @@ import .Microwave: radialcut!, radialstub! export Microwave, - bridge!, checkerboard!, device_template!, flux_bias!, diff --git a/src/microwave.jl b/src/microwave.jl index 4d3b7d6..22f37f2 100644 --- a/src/microwave.jl +++ b/src/microwave.jl @@ -9,7 +9,7 @@ using Devices.Rectangles using Devices.Polygons using Devices.Points using Devices.Cells -import Devices: AbstractPolygon, Coordinate, GDSMeta, Meta +import Devices: AbstractPolygon, Coordinate, GDSMeta, Meta, bridge! import Unitful: NoUnits, cm, μm, nm, ustrip, unit export bridge! diff --git a/src/paths/contstyles/compound.jl b/src/paths/contstyles/compound.jl index f033009..5ba6800 100644 --- a/src/paths/contstyles/compound.jl +++ b/src/paths/contstyles/compound.jl @@ -1,17 +1,18 @@ """ - struct CompoundStyle <: ContinuousStyle + struct CompoundStyle{T<:FloatCoordinate} <: ContinuousStyle{false} styles::Vector{Style} - grid::Vector{Float64} + grid::Vector{T} end Combines styles together, typically for use with a [`CompoundSegment`](@ref). -- `styles`: Array of styles making up the object. This is shallow-copied -by the outer constructor. +- `styles`: Array of styles making up the object. This is deep-copied by the outer + constructor. - `grid`: An array of `t` values needed for rendering the parameteric path. """ struct CompoundStyle{T<:FloatCoordinate} <: ContinuousStyle{false} styles::Vector{Style} grid::Vector{T} + tag::Symbol end function (s::CompoundStyle)(t) l0 = s.grid[1] @@ -23,16 +24,16 @@ function (s::CompoundStyle)(t) end return s.styles[length(s.grid) - 1], t - l0 end -copy(s::CompoundStyle) = (typeof(s))(deepcopy(s.styles), copy(s.grid)) -CompoundStyle(seg::AbstractVector, sty::AbstractVector) = - CompoundStyle(deepcopy(Vector{Style}(sty)), makegrid(seg, sty)) +copy(s::CompoundStyle, tag=gensym()) = (typeof(s))(deepcopy(s.styles), copy(s.grid), tag) +CompoundStyle(seg::AbstractVector{Segment{T}}, sty::AbstractVector, tag=gensym()) where {T} = + CompoundStyle(deepcopy(Vector{Style}(sty)), makegrid(seg, sty), tag) """ - makegrid{T<:Segment}(segments::AbstractVector{T}, styles) + makegrid(segments::AbstractVector{T}, styles) where {T<:Segment} Returns a collection with the values of `t` to use for rendering a `CompoundSegment` with a `CompoundStyle`. """ -function makegrid(segments::AbstractVector{T}, styles) where T<:Segment +function makegrid(segments::AbstractVector{T}, styles) where {T<:Segment} isempty(segments) && error("Cannot use makegrid with zero segments.") length(segments) != length(styles) && error("Must have same number of segments and styles.") @@ -52,3 +53,16 @@ for x in (:extent, :width) end summary(::CompoundStyle) = "Compound style" + +function translate(s::CompoundStyle, x, tag=gensym()) + s′ = copy(s, tag) + s′.grid .-= x + return s′ +end + +function pin(s::CompoundStyle; start=nothing, stop=nothing, tag=gensym()) + if start !== nothing + return translate(s, start, tag) + end + return copy(s, tag) +end diff --git a/src/paths/contstyles/cpw.jl b/src/paths/contstyles/cpw.jl index cc78697..dbeab5d 100644 --- a/src/paths/contstyles/cpw.jl +++ b/src/paths/contstyles/cpw.jl @@ -16,6 +16,7 @@ copy(x::GeneralCPW) = GeneralCPW(x.trace, x.gap) extent(s::GeneralCPW, t) = s.trace(t)/2 + s.gap(t) trace(s::GeneralCPW, t) = s.trace(t) gap(s::GeneralCPW, t) = s.gap(t) +translate(s::GeneralCPW, t) = GeneralCPW(x->s.trace(x+t), x->s.gap(x+t)) """ struct SimpleCPW{T<:Coordinate} <: CPW @@ -32,6 +33,7 @@ copy(x::SimpleCPW) = SimpleCPW(x.trace, x.gap) extent(s::SimpleCPW, t...) = s.trace/2 + s.gap trace(s::SimpleCPW, t...) = s.trace gap(s::SimpleCPW, t...) = s.gap +translate(s::SimpleCPW, t) = copy(s) """ CPW(trace::Coordinate, gap::Coordinate) diff --git a/src/paths/contstyles/decorated.jl b/src/paths/contstyles/decorated.jl index 4c77824..7a19df1 100644 --- a/src/paths/contstyles/decorated.jl +++ b/src/paths/contstyles/decorated.jl @@ -1,17 +1,17 @@ """ mutable struct DecoratedStyle{T<:FloatCoordinate} <: ContinuousStyle{false} s::Style - ts::Array{Float64,1} - dirs::Array{Int,1} - refs::Array{CellReference,1} + ts::Vector{Float64} + dirs::Vector{Int} + refs::Vector{CellReference} end Style with decorations, like structures periodically repeated along the path, etc. """ mutable struct DecoratedStyle{T<:FloatCoordinate} <: ContinuousStyle{false} s::Style - ts::Array{T,1} - dirs::Array{Int,1} - refs::Array{CellReference,1} + ts::Vector{T} + dirs::Vector{Int} + refs::Vector{CellReference} end summary(s::DecoratedStyle) = string(summary(s.s), " with ", length(s.refs), " decorations") @@ -85,3 +85,28 @@ function decorate(sty::DecoratedStyle, T, t, location, c) push!(sty.refs, c) sty end + +function pin(sty::DecoratedStyle{T}; start=nothing, stop=nothing) where T + x0 = ifelse(start === nothing, zero(sty.length), start) + x1 = ifelse(stop === nothing, sty.length, stop) + s = pin(undecorated(sty); start=start, stop=stop) + inds = findall(t->t in x0..x1, sty.ts) + ts = sty.ts[inds] + dirs = sty.dirs[inds] + refs = sty.refs[inds] + return DecoratedStyle{T}(s, ts, dirs, refs) +end + +""" + undecorate!(sty, t) +Removes all attachments at position `t` from a style. +""" +function undecorate!(sty::DecoratedStyle, t) + inds = findall(x->x==t, sty.ts) + deleteat!(sty.ts, inds) + deleteat!(sty.dirs, inds) + deleteat!(sty.refs, inds) + return sty +end +undecorate!(sty::Style, t) = sty +# handling compound style is probably brittle if it has decorations inside diff --git a/src/paths/contstyles/interface.jl b/src/paths/contstyles/interface.jl new file mode 100644 index 0000000..151204f --- /dev/null +++ b/src/paths/contstyles/interface.jl @@ -0,0 +1,35 @@ +""" + extent(s::Style, t) +For a style `s`, returns a distance tangential to the path specifying the lateral extent +of the polygons rendered. The extent is measured from the center of the path to the edge +of the polygon (half the total width along the path). The extent is evaluated at path length +`t` from the start of the associated segment. +""" +function extent end + +""" + width(s::Style, t) +For a style `s` and parameteric argument `t`, returns the width of paths rendered. +""" +function width end + +""" + translate(s::ContinuousStyle, x) +Creates a style `s′` such that all properties `f(s′, t) == f(s, t+x)`. Basically, advance +the style forward by path length `x`. +""" +function translate end + +""" + pin(s::ContinuousStyle; start=nothing, stop=nothing) +Imagine having a styled segment of length `L` split into two, the first segment having +length `l` and the second having length `L-l`. In all but the simplest styles, the styles +need to be modified in order to maintain the rendered appearances. A style appropriate for +the segment of length `l` (`L-l`) is given by `pin(s; stop=l)` (`pin(s; start=l)`). +""" +function pin(s::ContinuousStyle{false}; start=nothing, stop=nothing) + if start !== nothing + return translate(s, start) + end + return s +end diff --git a/src/paths/contstyles/strands.jl b/src/paths/contstyles/strands.jl index 407dcc9..72f6348 100644 --- a/src/paths/contstyles/strands.jl +++ b/src/paths/contstyles/strands.jl @@ -26,11 +26,14 @@ struct GeneralStrands{S,T,U} <: Strands{false} num::Int end copy(x::GeneralStrands) = GeneralStrands(x.offset, x.width, x.spacing, x.num) +GeneralStrands(o,w,s,n) = GeneralStrands{typeof(o), typeof(w), typeof(s)}(o,w,s,n) extent(s::GeneralStrands, t) = s.offset(t) + (s.num)*(s.width(t)) + (s.num - 1)*(s.spacing(t)) offset(s::GeneralStrands, t) = s.offset(t) width(s::GeneralStrands, t) = s.width(t) spacing(s::GeneralStrands, t) = s.spacing(t) num(s::GeneralStrands, t) = s.num +translate(s::GeneralStrands, t) = + GeneralStrands(x->s.offset(x+t), x->s.width(x+t), x->s.spacing(x+t), s.num) """ struct SimpleStrands{T<:Coordinate} <: Strands{false} @@ -62,6 +65,7 @@ offset(s::SimpleStrands, t...) = s.offset width(s::SimpleStrands, t...) = s.width spacing(s::SimpleStrands, t...) = s.spacing num(s::SimpleStrands, t...) = s.num +translate(s::SimpleStrands, t) = copy(s) """ Strands(offset::Coordinate, width::Coordinate, spacing::Coordinate, num::Int) diff --git a/src/paths/contstyles/tapers.jl b/src/paths/contstyles/tapers.jl index 3de4dc6..8f5efe2 100644 --- a/src/paths/contstyles/tapers.jl +++ b/src/paths/contstyles/tapers.jl @@ -16,6 +16,15 @@ end copy(x::TaperTrace) = TaperTrace(x.width_start, x_width_end) extent(s::TaperTrace, t) = 0.5 * width(s,t) width(s::TaperTrace, t) = (1-t/s.length) * s.width_start + t/s.length * s.width_end +function pin(s::TaperTrace; start=nothing, stop=nothing) + !isdefined(s, :length) && error("cannot `pin`; length of $s not yet determined.") + x0 = ifelse(start === nothing, zero(s.length), start) + x1 = ifelse(stop === nothing, s.length, stop) + @assert x1 - x0 > zero(x1 - x0) + @assert zero(x0) <= x0 < s.length + @assert zero(x1) < x1 <= s.length + return typeof(s)(width(s, x0), width(s, x1), x1-x0) +end function TaperTrace(width_start::Coordinate, width_end::Coordinate) dimension(width_start) != dimension(width_end) && throw(DimensionError(trace,gap)) w_s,w_e = promote(float(width_start), float(width_end)) @@ -55,6 +64,15 @@ function TaperCPW(trace_start::Coordinate, gap_start::Coordinate, trace_end::Coo t_s,g_s,t_e,g_e = promote(float(trace_start), float(gap_start), float(trace_end), float(gap_end)) return TaperCPW{typeof(t_s)}(t_s, g_s, t_e, g_e) end +function pin(sty::TaperCPW; start=nothing, stop=nothing) + !isdefined(sty, :length) && error("cannot `pin`; length of $sty not yet determined.") + x0 = ifelse(start === nothing, zero(sty.length), start) + x1 = ifelse(stop === nothing, sty.length, stop) + @assert x1 - x0 > zero(x1 - x0) + @assert zero(x0) <= x0 < sty.length + @assert zero(x1) < x1 <= sty.length + return typeof(sty)(trace(sty, x0), gap(sty, x0), trace(sty, x1), gap(sty, x1), x1 - x0) +end summary(s::TaperTrace) = string("Tapered trace with initial width ", s.width_start, " and final width ", s.width_end) diff --git a/src/paths/contstyles/trace.jl b/src/paths/contstyles/trace.jl index b19ccdb..400cc5b 100644 --- a/src/paths/contstyles/trace.jl +++ b/src/paths/contstyles/trace.jl @@ -12,6 +12,8 @@ end copy(x::GeneralTrace) = GeneralTrace(x.width) extent(s::GeneralTrace, t) = 0.5*s.width(t) width(s::GeneralTrace, t) = s.width(t) +trace(s::GeneralTrace, t) = s.width(t) +translate(s::GeneralTrace, t) = GeneralTrace(x->s.width(x+t)) """ struct SimpleTrace{T<:Coordinate} <: Trace{false} @@ -26,6 +28,8 @@ SimpleTrace(width) = SimpleTrace(width) copy(x::SimpleTrace) = Trace(x.width) extent(s::SimpleTrace, t...) = 0.5*s.width width(s::SimpleTrace, t...) = s.width +trace(s::SimpleTrace, t...) = s.width +translate(s::SimpleTrace, t) = copy(s) """ Trace(width) diff --git a/src/paths/intersect.jl b/src/paths/intersect.jl new file mode 100644 index 0000000..3905e96 --- /dev/null +++ b/src/paths/intersect.jl @@ -0,0 +1,297 @@ +using LinearAlgebra +import DataStructures: OrderedDict + +""" + intersectingpairs!(path1::Path, path2::Path) +Returns two vectors that can index into `path1` and `path2` respectively to return the +intersecting path nodes; also returns the intersection points. +""" +function intersectingpairs!(path1::Path, path2::Path) + T = promote_type(eltype(path1), eltype(path2)) + + # Approximate the paths by line segments and identify the intersecting line segments. + # If there are no intersections, we can return early. + inds1, segs1 = segmentize(path1) + local inds2, segs2, I, J + if path1 == path2 + # Want to ignore neighboring segments that intersect at their endpoints, + # so we'll use the one argument method of `intersectingpairs`. + inds2, segs2 = inds1, segs1 + I, J = intersectingpairs(segs1) + else + inds2, segs2 = segmentize(path2) + I, J = intersectingpairs(segs1, segs2) + end + # if `J` is empty, so is `I`, and there are no intersections. + isempty(J) && return inds1[I], inds2[J], Vector{Point{T}}() + + # `I`, `J` are indices into both `segs` and `inds` (respectively the + # `Polygons.LineSegment`s approximating the paths, and the indices of the + # `Paths.Segment`s they came from). Indices in `J` appear in ascending order and + # `I[j] < J[j]` for all valid `j`. + + i = 1 + intersections = Vector{Point{T}}(undef, length(J)) + for (s1, s2) in zip(segs1[I], segs2[J]) + # solve for the intersection of two lines. + l1, l2 = Polygons.Line(s1), Polygons.Line(s2) + intersections[i] = Polygons.intersection(l1, l2, false)[2] + i += 1 + end + + return inds1[I], inds2[J], intersections +end + +""" + segmentize(path::Paths.Path) +Turns a [`Paths.Path`](@ref) into an array of [`Polygons.LineSegment`](@ref) approximating +the path. Returns indices to mark which `LineSegment`s came from which [`Paths.Segment`](@ref). +""" +function segmentize(path::Paths.Path{T}) where T + segs0 = Polygons.LineSegment{T}[] + inds0 = Int[] + for (i,n) in enumerate(Paths.nodes(path)) + segs = segmentize(segment(n)) + inds = repeat([i], length(segs)) + append!(segs0, segs) + append!(inds0, inds) + end + return inds0, segs0 +end + +""" + segmentize(seg::Paths.Straight) +Returns a vector with a [`Polygons.LineSegment`](@ref) object corresponding to `seg`. +""" +function segmentize(seg::Paths.Straight) + return [Polygons.LineSegment(p0(seg), p1(seg))] +end + +function segmentize(seg::Paths.Corner{T}) where T + return Polygons.LineSegment{T}[] +end + +""" + segmentize(seg::Paths.Segment) +Generic fallback, approximating a [`Paths.Segment`](@ref) using many +[`Polygons.LineSegment`](@ref) objects. Returns a vector of `LineSegment`s. +""" +function segmentize(seg::Paths.Segment) + len = pathlength(seg) + bnds = (zero(len), len) + return segmentize(seg.(Devices.adapted_grid(t->Paths.direction(seg, t), bnds)), false) +end + +""" + intersectingpairs(segs1::Vector{<:Polygons.LineSegment}) +Returns two arrays `I`, `J` that index into `segs1` to return pairs of intersecting line +segments. It is guaranteed that `I[i] < J[i]` for all valid `i`, and that `J` is sorted +ascending. If any segments intersect at more than one point, an error is thrown. +""" +function intersectingpairs(segs1::Vector{<:Polygons.LineSegment}) + inds1 = Int[] + inds2 = Int[] + + s1 = collect(enumerate(segs1)) + for x in eachindex(s1) + (i,seg1) = s1[x] + for (j,seg2) in s1[(x+1):end] + intersects, atapoint, atendpoints = Polygons.intersects_at_endpoint(seg1, seg2) + intersects && !atapoint && error("degenerate segments: $seg1 and $seg2.") + if intersects && !atendpoints + # TODO: conceivable that we might want to handle atendpoints if the + # LineSegments came from non-adjacent Paths.Segments... + push!(inds1, i) + push!(inds2, j) + end + end + end + sortedinds = sortperm(inds2) + return inds1[sortedinds], inds2[sortedinds] +end + +""" + intersectingpairs(segs1::Vector{<:Polygons.LineSegment}, + segs2::Vector{<:Polygons.LineSegment}) +Returns two arrays `I`, `J` that index into `segs1` and `segs2` respectively to return +pairs of intersecting line segments. If any segments intersect at more than one point, an +error is thrown. It is guaranteed that `J` is sorted ascending. +""" +function intersectingpairs(segs1::Vector{<:Polygons.LineSegment}, + segs2::Vector{<:Polygons.LineSegment}) + inds1 = Int[] + inds2 = Int[] + + s1 = collect(enumerate(segs1)) + s2 = collect(enumerate(segs2)) + for (i,seg1) in s1 + for (j,seg2) in s2 + intersects, atapoint = Polygons.intersects(seg1, seg2) + intersects && !atapoint && error("degenerate segments: $seg1 and $seg2.") + if intersects + push!(inds1, i) + push!(inds2, j) + end + end + end + sortedinds = sortperm(inds2) + return inds1[sortedinds], inds2[sortedinds] +end + +""" + IntersectStyle{N} +Abstract type specifying "N-body interactions" for path intersection; `N::Int` typically 2. +""" +abstract type IntersectStyle{N} end + +""" + Bridges(; gap, footlength, bridgemetas, + bridgeprefix="bridges", bridgeunit=[nm or NoUnits]) +Style for automatically leaping one path over another with air bridges. + +The only required parameters to provide are: + +- `gap`: how much buffer to leave between the intersecting paths +- `footlength`: how long in the direction of the path is the bridge landing pad +- `bridgemetas`: a vector of `Meta` objects. The first element is for the bridge foot, + subsequent elements are for increasingly higher steps of the discretized bridge. The + length of this array therefore sets the vertical resolution. + +For each intersection, a bridge is rendered into a new cell. The cell has database units +`bridgeunit` which is `nm` by default (if `gap` has length units, otherwise no units). +Each cell is named uniquely with a prefix given by `bridgeprefix`. +""" +struct Bridges{S<:Meta, T<:Coordinate, U<:CoordinateUnits} <: IntersectStyle{2} + bridgemetas::Vector{S} + gap::T + footlength::T + bridgeprefix::String + bridgeunit::U + function Bridges(; gap::Coordinate, footlength::Coordinate, + bridgemetas::AbstractVector{<:Meta}, + bridgeprefix::AbstractString = "bridges", + bridgeunit::CoordinateUnits = (gap isa Length ? Unitful.nm : NoUnits)) + (dimension(gap) != dimension(footlength)) && + throw(DimensionError(gap, footlength)) + (dimension(gap) != dimension(bridgeunit)) && + throw(DimensionError(unit(gap), bridgeunit)) + + g, f = promote(gap, footlength) + S = eltype(bridgemetas) + T = eltype(g) + U = typeof(bridgeunit) + return new{S,T,U}(bridgemetas, g, f, bridgeprefix, bridgeunit) + end +end + +""" + intersect!(sty::IntersectStyle{2}, paths::Path...; + interactions::AbstractMatrix{Bool}=[default], reconcile=true) +Automatically modify paths to handle cases where they intersect. + +`interactions` is a keyword argument that specifies a pair-wise interaction matrix. If an +entry `(i,j)` is `true`, then `intersect_pairwise!(sty, paths[i], paths[j])` will be called. +By default the matrix is an upper-triangular matrix filled with `true`. This is typically +what you want: self-intersections of each path are handled, and the intersection of two +given paths is only handled once. + +Paths later in the argument list cross over paths earlier in the argument list. For +self-intersection (path with itself), segments later in a path will cross over segments +earlier in the same path (perhaps later this will be configurable by an option). +""" +function intersect!(sty::IntersectStyle{2}, paths::Path...; + interactions::AbstractMatrix{Bool} = + UpperTriangular(fill(true, (length(paths), length(paths)))), reconcile=true) + for I in findall(interactions) + intersect_pairwise!(sty, paths[I[1]], paths[I[2]]; reconcile=reconcile) + end +end + +""" + intersect_pairwise!(sty::Bridges, pa1::Path, pa2::Path; reconcile=true) +Automatically modify `pa2` to cross over `pa1` using air bridges. `pa1 === pa2` is acceptable. +The crossing-over segments must be `Paths.Straight`. +""" +function intersect_pairwise!(sty::Bridges{S}, pa1::Path, pa2::Path; reconcile=true) where S + # if `pa1 === pa2`, we don't want to disturb the indexing of `pa1` when we start + # splicing new segments into `pa2`... + _pa1 = (pa1 === pa2) ? copy(nodes(pa1)) : nodes(pa1) + + I, J, intersections = Paths.intersectingpairs!(pa1, pa2) + uJ = unique(J) + + # Bail immediately if any of the crossing-over segments are not straight + any(x->!(segment(x) isa Paths.Straight), view(pa2, uJ)) && + error("can only handle crossing of Straight paths.") + + # Compute the distances from the start of each segment identified in `J` to the + # corresponding intersection point(s). + dist2 = [norm(p - segment(pa2[j]).p0) for (j, p) in zip(J, intersections)] + + adj = 0 + for j in uJ # note that `uJ` is ascending because `J` was ascending. + # Each `j` in `uJ` is an index into `pa2` of a segment that needs to hop over another. + # The plan is to construct a new path consisting of only Paths.Straight and bridges + # and splice it into `pa2` at index `j`. To construct this new path, we need to + # start at the beginning of the existing straight segment at index `j` and + # identify which intersection points come first. We'll need the associated + # crossed-over segment indices too. + inds = findall(x->x==j, J) + increasing_order = sortperm(dist2[inds]) + Iinc, interinc = I[inds][increasing_order], intersections[inds][increasing_order] + + # We may have already spliced into the path, which will shift the indices from what + # they originally were. We keep track of this using `adj`. + j += adj + + # Get the old segment and style, and prepare a new path. + s2, st2 = segment(pa2[j]), style(pa2[j]) + pa = Path(p0(s2), α0=α0(s2)) + for (i, p, d) in zip(Iinc, interinc, dist2[inds][increasing_order]) + s1, st1 = segment(_pa1[i]), style(_pa1[i]) + + # Get the lateral extent of `pa1` at the intersection point in question. + # We should go straight up to the intersection point, less that extent, less + # the gap specified in `sty`. + howfar1 = norm(p - s1.p0) + extents = Paths.extent(st1, howfar1) + _howfar2 = d - extents - sty.gap + + # actually, we need to go a little bit less if we have a CPW, because we need + # to terminate the CPW... + termlA = Paths.terminationlength(st2, _howfar2) + howfar2 = _howfar2 - termlA + # TODO check that the straight section will be long enough to do this. + + # prepare a bridge of the correct size. + brc = Cell(uniquename(sty.bridgeprefix), sty.bridgeunit) + bridge!(brc, length(sty.bridgemetas) - 1, sty.footlength, + Paths.trace(st2, howfar2), 2*(sty.gap + extents + termlA), sty.bridgemetas) + + # now start going... + straight!(pa, howfar2 - pathlength(pa), st2) + (termlA > zero(termlA)) && + straight!(pa, termlA, Paths.terminationstyle(st2, _howfar2)) + nr = Paths.NoRender(2*Paths.extent(style1(pa), termlA)) + straight!(pa, sty.gap, nr) + straight!(pa, 2*extents, nr) + attach!(pa, CellReference(brc), extents) + straight!(pa, sty.gap, nr) + termlB = Paths.terminationlength(st2, pathlength(pa)) + (termlB > zero(termlB)) && + straight!(pa, termlB, Paths.terminationstyle(st2, pathlength(pa))) + end + # Finally, go the remaining length of the original segment. + straight!(pa, s2.l - pathlength(pa), st2) + + # Splice this path into `pa2`. Bam, automatic bridges. + splice!(pa2, j, pa; reconcile=false) + + # We've added this many segments to `pa2`, so adjust indices in the next pass... + adj += length(pa) - 1 + end + + # Finally, reconcile the path + # (adjust internal details for consistency beginning with the first spliced node). + reconcile && isempty(J) || reconcile!(pa2, first(J)) +end diff --git a/src/paths/norender.jl b/src/paths/norender.jl index 41eb493..a32a7f2 100644 --- a/src/paths/norender.jl +++ b/src/paths/norender.jl @@ -28,3 +28,6 @@ copy(x::SimpleNoRender) = SimpleNoRender(x.width) convert(::Type{DiscreteStyle}, x::NoRender) = NoRenderDiscrete() convert(::Type{ContinuousStyle}, x::NoRender) = NoRenderContinuous() NoRender(width::Coordinate) = SimpleNoRender(float(width)) + +translate(s::SimpleNoRender, t) = copy(s) +translate(s::NoRenderContinuous, t) = copy(s) diff --git a/src/paths/paths.jl b/src/paths/paths.jl index 690a434..0e1408d 100644 --- a/src/paths/paths.jl +++ b/src/paths/paths.jl @@ -14,6 +14,7 @@ import Base: empty!, deleteat!, length, + firstindex, lastindex, size, getindex, @@ -24,6 +25,9 @@ import Base: popfirst!, insert!, append!, + splice!, + split, + intersect!, show, summary, dims2string @@ -31,9 +35,11 @@ import Base: import Base.Iterators using ForwardDiff +import IntervalSets.(..) import Devices -import Devices: Coordinate, FloatCoordinate, GDSMeta, Meta -import Devices: bounds +import Devices: Polygons, Coordinate, FloatCoordinate, CoordinateUnits, GDSMeta, Meta +import Devices.Polygons: segmentize, intersects +import Devices: bounds, bridge! export Path @@ -56,24 +62,10 @@ export reconcile!, straight!, style, setstyle!, + terminate!, turn!, undecorated -""" - extent(s,t) -For a style `s`, returns a distance tangential to the path specifying the lateral extent -of the polygons rendered. The extent is measured from the center of the path to the edge -of the polygon (half the total width along the path). The extent is evaluated at path length -`t` from the start of the associated segment. -""" -function extent end - -""" - width(s,t) -For a style `s` and parameteric argument `t`, returns the width of paths rendered. -""" -function width end - """ abstract type Style end How to render a given path segment. @@ -97,6 +89,8 @@ Any style that applies to segments which have zero path length. """ abstract type DiscreteStyle <: Style end +include("contstyles/interface.jl") + """ abstract type Segment{T<:Coordinate} end Path segment in the plane. All Segment objects should have the implement the following @@ -301,17 +295,23 @@ end """ mutable struct Path{T<:Coordinate} <: AbstractVector{Node{T}} - p0::Point{T} - α0::Float64 - nodes::Array{Node{T},1} - end Type for abstracting an arbitrary styled path in the plane. Iterating returns -[`Paths.Node`](@ref) objects, essentially +[`Paths.Node`](@ref) objects. + + Path(p0::Point=Point(0.0,0.0); α0=0.0) + Path(p0x::Real, p0y::Real; kwargs...) + Path(p0::Point{T}; α0=0.0) where {T<:Length} + Path(p0x::T, p0y::T; kwargs...) where {T<:Length} + Path(p0x::Length, p0y::Length; kwargs...) + Path(u::LengthUnits; α0=0.0) + Path(v::Vector{<:Node}) +Convenience constructors for `Path{T}` object. """ mutable struct Path{T<:Coordinate} <: AbstractVector{Node{T}} p0::Point{T} α0::Float64 - nodes::Array{Node{T},1} + nodes::Vector{Node{T}} + laststyle::ContinuousStyle Path{T}() where {T} = new{T}(Point(zero(T),zero(T)), 0.0, Node{T}[]) Path{T}(a,b,c) where {T} = new{T}(a,b,c) @@ -333,20 +333,6 @@ function show(io::IO, x::Node) print(io, "$(segment(x)) styled as $(style(x))") end -""" -``` -Path(p0::Point=Point(0.0,0.0); α0=0.0) -Path(p0x::Real, p0y::Real; kwargs...) - -Path(p0::Point{T}; α0=0.0) where {T<:Length} -Path(p0x::T, p0y::T; kwargs...) where {T<:Length} -Path(p0x::Length, p0y::Length; kwargs...) - -Path(u::LengthUnits; α0=0.0) -``` - -Convenience constructors for `Path{T}` object. -""" function Path(p0::Point=Point(0.0,0.0); α0=0.0) Path{Float64}(p0, α0, Node{Float64}[]) end @@ -362,6 +348,11 @@ Path(p0x::Length, p0y::Length; kwargs...) = Path(promote(p0x,p0y)...; kwargs...) function Path(u::LengthUnits; α0=0.0) Path{typeof(0.0u)}(Point(0.0u,0.0u), α0, Node{typeof(0.0u)}[]) end +function Path(v::Vector{Node{T}}) where {T} + isempty(v) && return Path{T}(Point(zero(T), zero(T)), 0.0, v) + return Path{T}(p0(segment(v[1])), α0(segment(v[1])), v) +end + Path(x::Coordinate, y::Coordinate; kwargs...) = throw(DimensionError(x,y)) @@ -458,6 +449,8 @@ include("segments/turn.jl") include("segments/corner.jl") include("segments/compound.jl") +include("intersect.jl") + """ discretestyle1(p::Path) Returns the last-used discrete style in the path. @@ -466,9 +459,12 @@ discretestyle1(p::Path) = style1(p, DiscreteStyle) """ contstyle1(p::Path) -Returns the last-used discrete style in the path. +Returns the last user-provided continuous style in the path. """ -contstyle1(p::Path) = style1(p, ContinuousStyle) +function contstyle1(p::Path) + isdefined(p, :laststyle) || error("path is empty, provide a style.") + return p.laststyle +end """ reconcilelinkedlist!(p::Path, m::Integer) @@ -542,6 +538,18 @@ function reconcile!(p::Path, n::Integer=1) end # Methods for Path as AbstractArray + +function splice!(p::Path, inds; reconcile=true) + n = splice!(nodes(p), inds) + reconcile && reconcile!(p, first(inds)) + return n +end +function splice!(p::Path, inds, p2::Path; reconcile=true) + n = splice!(nodes(p), inds, nodes(p2)) + reconcile && reconcile!(p, first(inds)) + return n +end + length(p::Path) = length(nodes(p)) iterate(p::Path, state...) = iterate(nodes(p), state...) enumerate(p::Path) = enumerate(nodes(p)) @@ -551,143 +559,134 @@ Iterators.drop(p::Path, n::Int) = Iterators.drop(nodes(p), n) Iterators.cycle(p::Path) = Iterators.cycle(nodes(p)) isempty(p::Path) = isempty(nodes(p)) empty!(p::Path) = empty!(nodes(p)) -function deleteat!(p::Path, inds) +function deleteat!(p::Path, inds; reconcile=true) deleteat!(nodes(p), inds) - reconcile!(p, first(inds)) + reconcile && reconcile!(p, first(inds)) end +firstindex(p::Path) = 1 lastindex(p::Path) = length(nodes(p)) size(p::Path) = size(nodes(p)) getindex(p::Path, i::Integer) = nodes(p)[i] -function setindex!(p::Path, v::Node, i::Integer) +function setindex!(p::Path, v::Node, i::Integer; reconcile=true) nodes(p)[i] = v - reconcile!(p, i) + reconcile && reconcile!(p, i) end - -function setindex!(p::Path, v::Segment, i::Integer) - setsegment!(nodes(p)[i],v) - reconcile!(p, i) +function setindex!(p::Path, v::Segment, i::Integer; reconcile=true) + setsegment!(nodes(p)[i], v; reconcile=reconcile) end - -function setindex!(p::Path, v::Style, i::Integer) - setstyle!(nodes(p)[i],v) - reconcile!(p, i) +function setindex!(p::Path, v::Style, i::Integer; reconcile=true) + setstyle!(nodes(p)[i], v; reconcile=reconcile) end - -function push!(p::Path, node::Node) +function push!(p::Path, node::Node; reconcile=true) push!(nodes(p), node) - reconcile!(p, length(p)) + reconcile && reconcile!(p, length(p)) end - -function pushfirst!(p::Path, node::Node) +function pushfirst!(p::Path, node::Node; reconcile=true) pushfirst!(nodes(p), node) - reconcile!(p) + reconcile && reconcile!(p) end for x in (:push!, :pushfirst!) - @eval function ($x)(p::Path, seg::Segment, sty::Style) - ($x)(p, Node(seg,sty)) + @eval function ($x)(p::Path, seg::Segment, sty::Style; reconcile=true) + ($x)(p, Node(seg, sty); reconcile=reconcile) end - @eval function ($x)(p::Path, segsty0::Node, segsty::Node...) - ($x)(p, segsty0) + @eval function ($x)(p::Path, segsty0::Node, segsty::Node...; reconcile=true) + ($x)(p, segsty0; reconcile=reconcile) for x in segsty - ($x)(p, x) + ($x)(p, x; reconcile=reconcile) end end end -function pop!(p::Path) +function pop!(p::Path; reconcile=true) x = pop!(nodes(p)) - reconcile!(p, length(p)) - x + reconcile && reconcile!(p, length(p)) + return x end -function popfirst!(p::Path) +function popfirst!(p::Path; reconcile=true) x = popfirst!(nodes(p)) - reconcile!(p) - x + reconcile && reconcile!(p) + return x end -function insert!(p::Path, i::Integer, segsty::Node) +function insert!(p::Path, i::Integer, segsty::Node; reconcile=true) insert!(nodes(p), i, segsty) - reconcile!(p, i) + reconcile && reconcile!(p, i) end -insert!(p::Path, i::Integer, seg::Segment, sty::Style) = - insert!(p, i, Node(seg,sty)) +insert!(p::Path, i::Integer, seg::Segment, sty::Style; reconcile=true) = + insert!(p, i, Node(seg,sty); reconcile=reconcile) -function insert!(p::Path, i::Integer, seg::Segment) +function insert!(p::Path, i::Integer, seg::Segment; reconcile=true) if i == 1 sty = style0(p) else sty = style(nodes(p)[i-1]) end - insert!(p, i, Node(seg,sty)) + insert!(p, i, Node(seg,sty); reconcile=reconcile) end -function insert!(p::Path, i::Integer, segsty0::Node, segsty::Node...) - insert!(p, i, segsty0) +function insert!(p::Path, i::Integer, segsty0::Node, segsty::Node...; reconcile=true) + insert!(p, i, segsty0; reconcile=reconcile) for x in segsty - insert!(p, i, x) + insert!(p, i, x; reconcile=reconcile) end end """ - append!(p::Path, p′::Path) + append!(p::Path, p′::Path; reconcile=true) Given paths `p` and `p′`, path `p′` is appended to path `p`. The p0 and initial angle of the first segment from path `p′` is modified to match the last point and last angle of path `p`. """ -function append!(p::Path, p′::Path) +function append!(p::Path, p′::Path; reconcile=true) isempty(p′) && return i = length(p) lp, la = p1(p), α1(p) append!(nodes(p), nodes(p′)) - reconcile!(p, i+1) + reconcile && reconcile!(p, i+1) nothing end """ - simplify(p::Path, inds::UnitRange=1:length(p)) + simplify(p::Path, inds::UnitRange=firstindex(p):lastindex(p)) At `inds`, segments of a path are turned into a `CompoundSegment` and styles of a path are turned into a `CompoundStyle`. The method returns a tuple, `(segment, style)`. -- Indexing the path becomes more sane when you can combine several path -segments into one logical element. A launcher would have several indices -in a path unless you could simplify it. -- You don't need to think hard about boundaries between straights and turns -when you want a continuous styling of a very long path. +- Indexing the path becomes more sane when you can combine several path segments into one + logical element. A launcher would have several indices in a path unless you could simplify + it. +- You don't need to think hard about boundaries between straights and turns when you want a + continuous styling of a very long path. """ -function simplify(p::Path, inds::UnitRange=1:length(p)) - cseg = CompoundSegment(nodes(p)[inds]) - csty = CompoundStyle(cseg.segments, map(style, nodes(p)[inds])) +function simplify(p::Path, inds::UnitRange=firstindex(p):lastindex(p)) + tag = gensym() + cseg = CompoundSegment(nodes(p)[inds], tag) + csty = CompoundStyle(cseg.segments, map(style, nodes(p)[inds]), tag) Node(cseg, csty) end """ - simplify!(p::Path, inds::UnitRange=1:length(p)) + simplify!(p::Path, inds::UnitRange=firstindex(p):lastindex(p)) In-place version of [`simplify`](@ref). """ -function simplify!(p::Path, inds::UnitRange=1:length(p)) +function simplify!(p::Path, inds::UnitRange=firstindex(p):lastindex(p)) x = simplify(p, inds) deleteat!(p, inds) insert!(p, inds[1], x) p end -# function split{T<:Real}(s::CompoundSegment{T}, points) # WIP -# segs = CompoundSegment{T}[] -# segs -# end - """ meander!(p::Path, len, straightlen, r, α) -Alternate between going straight with length `straightlen` and turning -with radius `r` and angle `α`. Each turn goes the opposite direction of the -previous. The total length is `len`. Useful for making resonators. +Alternate between going straight with length `straightlen` and turning with radius `r` and +angle `α`. Each turn goes the opposite direction of the previous. The total length is `len`. +Useful for making resonators. -The straight and turn segments are combined into a `CompoundSegment` and -appended to the path `p`. +The straight and turn segments are combined into a `CompoundSegment` and appended to the +path `p`. """ function meander!(p::Path, len, straightlen, r, α) unit = straightlen + r*abs(α) @@ -753,4 +752,86 @@ function _launch!(p::Path{T}; kwargs...) where {T <: Coordinate} CPW(trace1, gap1) end +""" + terminate!(pa::Path) +End a path with open termination (do nothing for a trace, leave a gap for a CPW). +""" +function terminate!(pa::Path{T}) where T + terminationlength(pa) > zero(T) && + straight!(pa, terminationlength(pa), terminationstyle(pa)) +end + +terminationlength(pa::Path, t=pathlength(pa[end])) = terminationlength(style(pa[end]), t) +terminationlength(s::Trace, t) = zero(t) +terminationlength(s::CPW, t) = gap(s,t) + +terminationstyle(pa::Path, t=pathlength(pa[end])) = terminationstyle(style(pa[end]), t) +terminationstyle(s::CPW, t) = Paths.Trace(2 * extent(s,t)) + +""" + split(n::Node, x) +Splits a path node at a position `x` along the segment, returning a path. + +Splitting and splicing back into a path: + splice!(path, i, split(path[i], x)) +""" +function split(n::Node, x::Coordinate) + seg1, seg2, sty1, sty2 = split(segment(n), style(n), x) + + n1 = Node(seg1, sty1) + n2 = Node(seg2, sty2) + n1.prev = n.prev + n1.next = n2 + n2.prev = n1 + n2.next = n.next + + return Path([n1, n2]) +end + +""" + split(n::Node, x::AbstractVector{<:Coordinate}; issorted=false) +Splits a path node at positions in `x` along the segment, returning a path. +If `issorted`, don't sort `x` first (otherwise required for this to work). +""" +function split(n::Node, x::AbstractVector{<:Coordinate}; issorted=false) + @assert !isempty(x) + sortedx = issorted ? x : sort(x) + + i = 2 + L = first(sortedx) + path = split(n, L) + for pos in view(sortedx, (firstindex(sortedx)+1):lastindex(sortedx)) + splice!(path, i, split(path[i], pos-L); reconcile=false) + L = pos + i += 1 + end + reconcile!(path) + return path +end + +function split(seg::Segment, sty::Style, x) + return (split(seg, x)..., split(sty, x)...) +end + +function split(seg::Segment, x) + @assert zero(x) < x < pathlength(seg) + @assert seg isa ContinuousSegment + return _split(seg, x) +end + +function split(sty::Style, x) + @assert sty isa ContinuousStyle + return _split(sty, x) +end + +function _split(sty::Style, x) + s1, s2 = pin(sty; stop=x), pin(sty; start=x) + undecorate!(s1, x) # don't duplicate attachments at the split point! + return s1, s2 +end + +function _split(sty::CompoundStyle, x, tag1, tag2) + return pin(sty; stop=x, tag=tag1), pin(sty; start=x, tag=tag2) +end + end diff --git a/src/paths/segments/compound.jl b/src/paths/segments/compound.jl index 77a7376..ae4cb4b 100644 --- a/src/paths/segments/compound.jl +++ b/src/paths/segments/compound.jl @@ -11,16 +11,22 @@ path function, and are not allowed in a `CompoundSegment`. """ struct CompoundSegment{T} <: ContinuousSegment{T} segments::Vector{Segment{T}} + tag::Symbol - CompoundSegment{T}(segments) where {T} = begin + function CompoundSegment{T}(segments, tag=gensym()) where {T} if any(x->isa(x,Corner), segments) error("cannot have corners in a `CompoundSegment`. You may have ", "tried to simplify a path containing `Corner` objects.") else - new{T}(deepcopy(Array(segments))) + new{T}(deepcopy(Array(segments)), tag) end end end +copy(s::CompoundSegment) = (typeof(s))(s.segments) +CompoundSegment(nodes::AbstractVector{Node{T}}, tag=gensym()) where {T} = + CompoundSegment{T}(map(segment, nodes), tag) +CompoundSegment(segments::AbstractVector{Segment{T}}, tag=gensym()) where {T} = + CompoundSegment{T}(segments, tag) # Parametric function over the domain [zero(T),pathlength(c)] that represents the # compound segments. @@ -58,11 +64,35 @@ function (s::CompoundSegment{T})(t) where {T} end end -CompoundSegment(nodes::AbstractVector{Node{T}}) where {T} = - CompoundSegment{T}(map(segment, nodes)) +function _split(seg::CompoundSegment{T}, x, tag1=gensym(), tag2=gensym()) where {T} + @assert zero(x) < x < pathlength(seg) + c = seg.segments + isempty(c) && error("cannot split a CompoundSegment based on zero segments.") + + L = pathlength(seg) + l0 = zero(L) + + for i in firstindex(c):lastindex(c) + seg = c[i] + l1 = l0 + pathlength(seg) + if l0 <= x < l1 + if x == l0 # can't happen on the firstindex because we have an assertion earlier + # This is a clean split between segments + seg1 = CompoundSegment(c[firstindex(c):(i-1)], tag1) + seg2 = CompoundSegment(c[i:lastindex(c)], tag2) + return seg1, seg2 + else + s1, s2 = split(seg, x - l0) + seg1 = CompoundSegment(push!(c[firstindex(c):(i-1)], s1), tag1) + seg2 = CompoundSegment(pushfirst!(c[(i+1):lastindex(c)], s2), tag2) + return seg1, seg2 + end + end + l0 = l1 + end +end summary(s::CompoundSegment) = string(length(s.segments), " segments") -copy(s::CompoundSegment) = (typeof(s))(s.segments) pathlength(s::CompoundSegment) = sum(pathlength, s.segments) function setα0p0!(s::CompoundSegment, angle, p::Point) diff --git a/src/paths/segments/corner.jl b/src/paths/segments/corner.jl index fb2b992..713d26d 100644 --- a/src/paths/segments/corner.jl +++ b/src/paths/segments/corner.jl @@ -26,13 +26,14 @@ copy(x::Corner{T}) where {T} = Corner{T}(x.α, x.p0, x.α0, x.extent) pathlength(::Corner{T}) where {T} = zero(T) p0(s::Corner) = s.p0 -function p1(s::Corner) - sgn = ifelse(s.α >= 0.0, 1, -1) - ∠A = s.α0+sgn*π/2 - v = s.extent*Point(cos(∠A),sin(∠A)) - ∠B = ∠A + π + s.α - s.p0+v+s.extent*Point(cos(∠B),sin(∠B)) -end +p1(s::Corner) = s.p0 +# function p1(s::Corner) +# sgn = ifelse(s.α >= 0.0, 1, -1) +# ∠A = s.α0+sgn*π/2 +# v = s.extent*Point(cos(∠A),sin(∠A)) +# ∠B = ∠A + π + s.α +# s.p0+v+s.extent*Point(cos(∠B),sin(∠B)) +# end α0(s::Corner) = s.α0 α1(s::Corner) = s.α0 + s.α setp0!(s::Corner, p::Point) = s.p0 = p @@ -47,8 +48,17 @@ The style chosen for this corner, if not specified, is the last `DiscreteStyle` path. """ function corner!(p::Path, α, sty::Style=discretestyle1(p)) + segment(last(p)) isa Paths.Straight || + error("corners must follow `Paths.Straight` segments.") T = eltype(p) seg = Corner{T}(α) + + ext = extent(style(p[end]), pathlength(p[end])) + w = ext * tan(abs(0.5*seg.α)) + pa = split(p[end], pathlength(p[end]) - w) + setstyle!(pa[end], SimpleNoRender(2*ext)) + splice!(p, length(p), pa) + # convert takes NoRender() → NoRenderDiscrete() push!(p, Node(seg, convert(DiscreteStyle, sty))) nothing diff --git a/src/paths/segments/straight.jl b/src/paths/segments/straight.jl index df6aca4..b534ef0 100644 --- a/src/paths/segments/straight.jl +++ b/src/paths/segments/straight.jl @@ -49,7 +49,28 @@ function straight!(p::Path, l::Coordinate, sty::Style=contstyle1(p)) dimension(T) != dimension(typeof(l)) && throw(DimensionError(T(1),l)) @assert l >= zero(l) "tried to go straight by a negative amount." seg = Straight{T}(l, p1(p), α1(p)) - # convert takes NoRender() → NoRenderContinuous() - push!(p, Node(seg, convert(ContinuousStyle, sty))) + + if !isempty(p) && (segment(last(p)) isa Paths.Corner) + cseg = segment(last(p)) + minlen = cseg.extent * tan(abs(0.5*cseg.α)) + @assert l > minlen "straight following corner needs minimum length $l." + + # convert takes NoRender() → NoRenderContinuous() + push!(p, Node(seg, convert(ContinuousStyle, sty))) + + pa = split(p[end], minlen) + setstyle!(pa[1], SimpleNoRender(2*cseg.extent)) + splice!(p, length(p), pa) + else + # convert takes NoRender() → NoRenderContinuous() + push!(p, Node(seg, convert(ContinuousStyle, sty))) + end + p.laststyle = sty nothing end + +function _split(seg::Straight{T}, x) where {T} + s1 = Straight{T}(x, seg.p0, seg.α0) + s2 = Straight{T}(seg.l - x, seg(x), seg.α0) + return s1, s2 +end diff --git a/src/paths/segments/turn.jl b/src/paths/segments/turn.jl index d056aff..082372f 100644 --- a/src/paths/segments/turn.jl +++ b/src/paths/segments/turn.jl @@ -10,9 +10,11 @@ mutable struct Turn{T} <: ContinuousSegment{T} α0::Float64 end function (s::Turn)(t) + # assuming s.α not zero x = ifelse(s.r == zero(s.r), typeof(s.r)(one(s.r)), s.r) # guard against div by zero - cen = s.p0 + Point(s.r*cos(s.α0+sign(s.α)*π/2), s.r*sin(s.α0+sign(s.α)*π/2)) - cen + Point(s.r*cos(s.α0-sign(s.α)*(π/2-t/x)), s.r*sin(s.α0-sign(s.α)*(π/2-t/x))) + cen = s.p0 + Point(-s.r*sign(s.α)sin(s.α0), s.r*sign(s.α)cos(s.α0)) + return cen + Point(s.r*(cos(s.α0)*sin(t/x) + sign(s.α)*sin(s.α0)*cos(t/x)), + s.r*(sin(s.α0)*sin(t/x) - sign(s.α)*cos(s.α0)*cos(t/x))) end """ @@ -52,8 +54,11 @@ Positive angle turns left. By default, we take the last continuous style in the function turn!(p::Path, α, r::Coordinate, sty::Style=contstyle1(p)) T = eltype(p) dimension(T) != dimension(typeof(r)) && throw(DimensionError(T(1),r)) + !isempty(p) && (segment(last(p)) isa Paths.Corner) && + error("`Paths.Straight` segments must follow `Paths.Corner`s.") seg = Turn{T}(α, r, p1(p), α1(p)) push!(p, Node(seg, convert(ContinuousStyle, sty))) + p.laststyle = sty nothing end @@ -70,6 +75,8 @@ By default, we take the last continuous style in the path. function turn!(p::Path, str::String, r::Coordinate, sty::Style=contstyle1(p)) T = eltype(p) dimension(T) != dimension(typeof(r)) && throw(DimensionError(T(1),r)) + !isempty(p) && (segment(last(p)) isa Paths.Corner) && + error("`Paths.Straight` segments must follow `Paths.Corner`s.") for ch in str if ch == 'l' α = π/2 @@ -82,5 +89,15 @@ function turn!(p::Path, str::String, r::Coordinate, sty::Style=contstyle1(p)) # convert takes NoRender() → NoRenderContinuous() push!(p, Node(seg, convert(ContinuousStyle, sty))) end + p.laststyle = sty nothing end + +function _split(seg::Turn{T}, x) where {T} + r, α = seg.r, seg.α + α1 = x / r + α2 = α - α1 + s1 = Turn{T}(α1, seg.r, seg.p0, seg.α0) + s2 = Turn{T}(α2, seg.r, seg(x), seg.α0 + α1) + return s1, s2 +end diff --git a/src/polygons.jl b/src/polygons.jl index 9cfee80..3bd117e 100644 --- a/src/polygons.jl +++ b/src/polygons.jl @@ -9,18 +9,21 @@ import Base: copy, promote_rule using ForwardDiff import CoordinateTransformations: AffineMap, LinearMap, Translation import Clipper -import Clipper: orientation, children, contour +import Clipper: children, contour import StaticArrays import Devices import Devices: AbstractPolygon, Coordinate, GDSMeta, Meta -import Devices: bounds, lowerleft, upperright +import Devices: bounds, lowerleft, upperright, orientation import Unitful import Unitful: Quantity, Length, dimension, unit, ustrip, uconvert, ° using ..Points using ..Rectangles import ..cclipper +import IntervalSets.(..) +import IntervalSets.endpoints + export Polygon export points export clip, offset, circle @@ -424,37 +427,70 @@ end ### cutting algorithm -abstract type D1 end +abstract type D1{T} end +Δy(d1::D1) = d1.p1.y - d1.p0.y +Δx(d1::D1) = d1.p1.x - d1.p0.x ab(p0, p1) = Point(gety(p1)-gety(p0), getx(p0)-getx(p1)) -struct Segment <: D1 - p0::Point{Float64} - p1::Point{Float64} - ab::Point{Float64} +""" + LineSegment{T} <: D1{T} +Represents a line segment. By construction, `p0.x <= p1.x`. +""" +struct LineSegment{T} <: D1{T} + p0::Point{T} + p1::Point{T} + function LineSegment(p0::Point{T}, p1::Point{T}) where T + if p1.x < p0.x + return new{T}(p1, p0) + else + return new{T}(p0, p1) + end + end end -Segment(p0,p1) = Segment(p0, p1, ab(p0, p1)) +LineSegment(p0::Point{S}, p1::Point{T}) where {S,T} = LineSegment(promote(p0, p1)...) -struct Ray <: D1 - p0::Point{Float64} - p1::Point{Float64} - ab::Point{Float64} +""" + Ray{T} <: D1{T} +Represents a ray. The ray starts at `p0` and goes toward `p1`. +""" +struct Ray{T} <: D1{T} + p0::Point{T} + p1::Point{T} end -Ray(p0,p1) = Ray(p0, p1, ab(p0, p1)) +Ray(p0::Point{S}, p1::Point{T}) where {S,T} = Ray(promote(p0, p1)...) -struct Line <: D1 - p0::Point{Float64} - p1::Point{Float64} - ab::Point{Float64} +struct Line{T} <: D1{T} + p0::Point{T} + p1::Point{T} end -Line(p0,p1) = Line(p0, p1, ab(p0, p1)) +Line(p0::Point{S}, p1::Point{T}) where {S,T} = Line(promote(p0, p1)...) +Line(seg::LineSegment) = Line(seg.p0, seg.p1) + +Base.promote_rule(::Type{Line{S}}, ::Type{Line{T}}) where {S,T} = Line{promote_type(S,T)} +Base.convert(::Type{Line{S}}, L::Line) where S = Line{S}(L.p0, L.p1) -function segments(vertices) +""" + segmentize(vertices, closed=true) +Make an array of `LineSegment` out of an array of points. If `closed`, a segment should go +between the first and last point, otherwise nah. +""" +function segmentize(vertices, closed=true) l = length(vertices) - [Segment(vertices[i], vertices[i==l ? 1 : i+1]) for i = 1:l] + if closed + return [LineSegment(vertices[i], vertices[i==l ? 1 : i+1]) for i = 1:l] + else + return [LineSegment(vertices[i], vertices[i+1]) for i = 1:(l-1)] + end end -# Find the lower-most then left-most polygon +""" + uniqueray(v::Vector{Point{T}}) where {T <: Real} +Given an array of points (thought to indicate a polygon or a hole in a polygon), +find the lowest / most negative y-coordinate[s] `miny`, then the lowest / most negative +x-coordinate `minx` of the points having that y-coordinate. This `Point(minx,miny)` ∈ `v`. +Return a ray pointing in -ŷ direction from that point. +""" function uniqueray(v::Vector{Point{T}}) where {T <: Real} nopts = reinterpret(T, v) yarr = view(nopts, 2:2:length(nopts)) @@ -464,42 +500,127 @@ function uniqueray(v::Vector{Point{T}}) where {T <: Real} Ray(Point(minx,miny), Point(minx, miny-1)) end -orientation(p::Polygon) = orientation(reinterpret(Clipper.IntPoint, p.p)) +""" + orientation(p::Polygon) +Returns 1 if the points in the polygon contour are going counter-clockwise, -1 if clockwise. +Clipper considers clockwise-oriented polygons to be holes for some polygon fill types. +""" +function orientation(p::Polygon) + ccall((:orientation, cclipper), Cuchar, (Ptr{Clipper.IntPoint}, Csize_t), + reinterpret(Clipper.IntPoint, p.p),length(p.p)) == 1 ? 1 : -1 +end + +""" + ishole(p::Polygon) +Returns `true` if Clipper would consider this polygon to be a hole, for applicable +polygon fill rules. +""" +ishole(p::Polygon) = orientation(p) == -1 + +""" + orientation(p1::Point, p2::Point, p3::Point) +Returns 1 if the path `p1`--`p2`--`p3` is going counter-clockwise (increasing angle), +-1 if the path is going clockwise (decreasing angle), 0 if `p1`, `p2`, `p3` are colinear. +""" +function orientation(p1::Point, p2::Point, p3::Point) + return sign((p3.y-p2.y)*(p2.x-p1.x)-(p2.y-p1.y)*(p3.x-p2.x)) +end -ishole(p) = orientation(p) == false -isparallel(A::D1, B::D1) = getx(A.ab) * gety(B.ab) == getx(B.ab) * gety(A.ab) -isdegenerate(A::D1, B::D1) = dot(A.ab, B.p0-A.p0) == dot(A.ab, B.p1-A.p0) == 0 +isparallel(A::D1, B::D1) = Δy(A) * Δx(B) == Δy(B) * Δx(A) +isdegenerate(A::D1, B::D1) = + orientation(A.p0, A.p1, B.p0) == orientation(A.p0, A.p1, B.p1) == 0 +iscolinear(A::D1, B::Point) = orientation(A.p0, A.p1, B) == orientation(B, A.p1, A.p0) == 0 +iscolinear(A::Point, B::D1) = iscolinear(B, A) -# Expected to be fast -function intersects(A::Segment, B::Segment) - sb0 = sign(dot(A.ab, B.p0-A.p0)) - sb1 = sign(dot(A.ab, B.p1-A.p0)) +""" + intersects(A::LineSegment, B::LineSegment) +Returns two `Bool`s: +1) Does `A` intersect `B`? +2) Did an intersection happen at a single point? (`false` if no intersection) +""" +function intersects(A::LineSegment, B::LineSegment) + sb0 = orientation(A.p0, A.p1, B.p0) + sb1 = orientation(A.p0, A.p1, B.p1) sb = sb0 == sb1 - sa0 = sign(dot(B.ab, A.p0-B.p0)) - sa1 = sign(dot(B.ab, A.p1-B.p0)) + sa0 = orientation(B.p0, B.p1, A.p0) + sa1 = orientation(B.p0, B.p1, A.p1) sa = sa0 == sa1 if sa == false && sb == false - return true + return true, true else - return false + # Test for special case of colinearity + if sb0 == sb1 == sa0 == sa1 == 0 + xinter = intersect(A.p0.x..A.p1.x, B.p0.x..B.p1.x) + yinter = intersect(A.p0.y..A.p1.y, B.p0.y..B.p1.y) + if !isempty(xinter) && !isempty(yinter) + if reduce(==, endpoints(xinter)) && reduce(==, endpoints(yinter)) + return true, true + else + return true, false + end + else + return false, false + end + else + return false, false + end + end +end + +""" + intersects_at_endpoint(A::LineSegment, B::LineSegment) +Returns three `Bool`s: +1) Does `A` intersect `B`? +2) Did an intersection happen at a single point? (`false` if no intersection) +3) Did an endpoint of `A` intersect an endpoint of `B`? +""" +function intersects_at_endpoint(A::LineSegment, B::LineSegment) + A_intersects_B, atapoint = intersects(A,B) + if A_intersects_B + if atapoint + if (A.p1 == B.p0) || (A.p1 == B.p1) || (A.p0 == B.p0) || (A.p0 == B.p1) + return A_intersects_B, atapoint, true + else + return A_intersects_B, atapoint, false + end + else + return A_intersects_B, atapoint, false + end + else + return A_intersects_B, atapoint, false end end -function onray(p::Point{T}, A::Ray) where {T <: Real} - return (dot(A.ab, p-A.p0) ≈ 0) && - (dot(A.p1-A.p0, p-A.p0) >= 0) +""" + intersects(p::Point, A::Ray) +Does `p` intersect `A`? +""" +function intersects(p::Point, A::Ray) + correctdir = dot(A.p1-A.p0, p-A.p0) >= 0 + return iscolinear(p, A) && correctdir end -function onsegment(p::Point{T}, A::Segment) where {T <: Real} - return (dot(A.ab, p-A.p0) ≈ 0) && - (dot(A.p1-A.p0, p-A.p0) >= 0) && - (dot(A.p0-A.p1, p-A.p1) >= 0) +""" + intersects(p::Point, A::LineSegment) +Does `p` intersect `A`? +""" +function intersects(p::Point, A::LineSegment) + if iscolinear(p, A) + xinter = intersect(A.p0.x..A.p1.x, p.x..p.x) + yinter = intersect(A.p0.y..A.p1.y, p.y..p.y) + if !isempty(xinter) && !isempty(yinter) + return true + else + return false + end + else + return false + end end -# Not type stable... -function intersection(A::Ray, B::Segment) +function intersection(A::Ray, B::LineSegment) if isparallel(A, B) if isdegenerate(A, B) # correct direction? @@ -525,31 +646,45 @@ function intersection(A::Ray, B::Segment) return false, Point(0.,0.) end else - tf, where = intersection(Line(A.p0,A.p1,A.ab), Line(B.p0,B.p1,B.ab), false) - if onray(where, A) && onsegment(where, B) - return true, where + tf, w = intersection(Line(A.p0,A.p1), Line(B.p0,B.p1), false) + if intersects(w, A) && intersects(w, B) + return true, w else return false, Point(0.,0.) end end end -function intersection(A::Line, B::Line, checkparallel=true) +function intersection(A::Line{T}, B::Line{T}, checkparallel=true) where T if checkparallel # parallel checking goes here! else - w = [getx(A.ab) gety(A.ab); getx(B.ab) gety(B.ab)] \ [dot(A.ab, A.p0), dot(B.ab, B.p0)] - true, Point(w) + xA, xB, yA, yB = Δx(A), Δx(B), Δy(A), Δy(B) + w = ustrip.([yA -xA; yB -xB]) \ ustrip.([A.p0.x*yA - A.p0.y*xA, B.p0.x*yB - B.p0.y*xB]) + true, Point(w)*unit(A.p0.x) end end +""" + interiorcuts(nodeortree::Clipper.PolyNode, outpolys::Vector{Polygon{T}}) where {T} +Clipper gives polygons with holes as separate contours. The GDS-II format doesn't support +this. This function makes cuts between the inner/outer contours so that ultimately there +is just one contour with one or more overlapping edges. + +Example: +┌────────────┐ ┌────────────┐ +│ ┌──┐ │ becomes... │ ┌──┐ │ +│ └──┘ ┌──┐ │ │ ├──┘ ┌──┐ │ +│ └──┘ │ │ │ ├──┘ │ +└────────────┘ └─┴─────┴────┘ +""" function interiorcuts(nodeortree::Clipper.PolyNode, outpolys::Vector{Polygon{T}}) where {T} # Assumes we have first element an enclosing polygon with the rest being holes. # We also assume no hole collision. minpt = Point(-Inf, -Inf) for enclosing in children(nodeortree) - segs = segments(contour(enclosing)) + segs = segmentize(contour(enclosing)) for hole in children(enclosing) # process all the holes. interiorcuts(hole, outpolys) @@ -561,29 +696,20 @@ function interiorcuts(nodeortree::Clipper.PolyNode, outpolys::Vector{Polygon{T}} k = -1 bestwhere = minpt for (j,s) in enumerate(segs) - tf, where = intersection(ray, s) + tf, wh = intersection(ray, s) if tf - if gety(where) > gety(bestwhere) - bestwhere = where + if gety(wh) > gety(bestwhere) + bestwhere = wh k = j end end end - # println(bestwhere) - # println(k) - # println(ray) # Since the polygon was enclosing, an intersection had to happen *somewhere*. if k != -1 w = Point{Int64}(round(getx(bestwhere)), round(gety(bestwhere))) - # println(w) kp1 = contour(enclosing)[(k+1 > length(contour(enclosing))) ? 1 : k+1] - # println(contour(enclosing)[1:k]) - # println(w) - # println(contour(hole)) - # println(contour(enclosing)[(k+1):end]) - # Make the cut in the enclosing polygon enclosing.contour = Point{Int64}[contour(enclosing)[1:k]; [w]; # could actually be contour(enclosing)[k]... should check for this. @@ -594,11 +720,11 @@ function interiorcuts(nodeortree::Clipper.PolyNode, outpolys::Vector{Polygon{T}} # update the segment cache segs = [segs[1:(k-1)]; - Segment(contour(enclosing)[k], w); - Segment(w, contour(hole)[1]); - segments(contour(hole)); - Segment(contour(hole)[1], w); - Segment(w, kp1); + LineSegment(contour(enclosing)[k], w); + LineSegment(w, contour(hole)[1]); + segmentize(contour(hole)); + LineSegment(contour(hole)[1], w); + LineSegment(w, kp1); segs[(k+1):end]] end end diff --git a/src/predicates.jl b/src/predicates.jl new file mode 100644 index 0000000..0fd57cc --- /dev/null +++ b/src/predicates.jl @@ -0,0 +1,275 @@ +#= +See LICENSE.md. +=# + +const Float64Like = Union{Float64, Length{Float64}} + +macro Fast_Two_Sum_Tail(a, b, x, y) + return esc(quote + bvirt = $x - $a + $y = $b - bvirt + end) +end + +macro Fast_Two_Sum(a, b, x, y) + return esc(quote + $x = $a + $b + @Fast_Two_Sum_Tail($a, $b, $x, $y) + end) +end + +macro Two_Sum_Tail(a, b, x, y) + return esc(quote + bvirt = $x - $a + avirt = $x - bvirt + bround = $b - bvirt + around = $a - avirt + $y = around + bround + end) +end + +macro Two_Sum(a, b, x, y) + return esc(quote + $x = $a + $b + @Two_Sum_Tail($a, $b, $x, $y) + end) +end + +macro Two_Diff_Tail(a, b, x, y) + return esc(quote + bvirt = $a - $x + avirt = $x + bvirt + bround = bvirt - $b + around = $a - avirt + $y = around + bround + end) +end + +macro Two_Diff(a, b, x, y) + return esc(quote + $x = $a - $b + @Two_Diff_Tail($a, $b, $x, $y) + end) +end + +macro Split(a, ahi, alo) + return esc(quote + c = splitter * $a + abig = c - $a + $ahi = c - abig + $alo = $a - $ahi + end) +end + +macro Two_One_Diff(a1, a0, b, x2, x1, x0) + return esc(quote + @Two_Diff($a0, $b , _i, $x0) + @Two_Sum( $a1, _i, $x2, $x1) + end) +end + +macro Two_Two_Diff(a1, a0, b1, b0, x3, x2, x1, x0) + return esc(quote + @Two_One_Diff($a1, $a0, $b0, _j, _0, $x0) + @Two_One_Diff(_j, _0, $b1, $x3, $x2, $x1) + end) +end + +macro Two_Product_Tail(a, b, x, y) + return esc(quote + @Split($a, ahi, alo) + @Split($b, bhi, blo) + err1 = $x - (ahi * bhi) + err2 = err1 - (alo * bhi) + err3 = err2 - (ahi * blo) + $y = (alo * blo) - err3 + end) +end + +macro Two_Product(a, b, x, y) + return esc(quote + $x = $a * $b + @Two_Product_Tail($a, $b, $x, $y) + end) +end + +function orient2dadapt(pa::Point{T}, pb::Point{T}, pc::Point{T}, detsum, + ccwerrboundB::Float64 = Devices.ccwerrboundB, + ccwerrboundC::Float64 = Devices.ccwerrboundC, + resulterrbound::Float64 = Devices.resulterrbound) where {T} + z = zero(T)*zero(T) + B = Vector{typeof(z)}(undef, 4) + C1 = Vector{typeof(z)}(undef, 8) + C2 = Vector{typeof(z)}(undef, 12) + D = Vector{typeof(z)}(undef, 16) + u = Vector{typeof(z)}(undef, 4) + + acx = pa[1] - pc[1] + bcx = pb[1] - pc[1] + acy = pa[2] - pc[2] + bcy = pb[2] - pc[2] + + @Two_Product(acx, bcy, detleft, detlefttail) + @Two_Product(acy, bcx, detright, detrighttail) + + @Two_Two_Diff(detleft, detlefttail, detright, detrighttail, + B3, B[3], B[2], B[1]) + B[4] = B3 + + det = sum(B) + errbound = ccwerrboundB * detsum + if (det >= errbound) || (-det >= errbound) + return sign(det) + end + + @Two_Diff_Tail(pa[1], pc[1], acx, acxtail) + @Two_Diff_Tail(pb[1], pc[1], bcx, bcxtail) + @Two_Diff_Tail(pa[2], pc[2], acy, acytail) + @Two_Diff_Tail(pb[2], pc[2], bcy, bcytail) + + if (acxtail == z) && (acytail == z) && (bcxtail == z) && (bcytail == z) + return sign(det) + end + + errbound = ccwerrboundC * detsum + resulterrbound * abs(det) + det += (acx * bcytail + bcy * acxtail) - (acy * bcxtail + bcx * acytail) + if (det >= errbound) || (-det >= errbound) + return sign(det) + end + + @Two_Product(acxtail, bcy, s1, s0) + @Two_Product(acytail, bcx, t1, t0) + @Two_Two_Diff(s1, s0, t1, t0, u3, u[3], u[2], u[1]) + u[4] = u3 + C1length = fast_expansion_sum_zeroelim(4, B, 4, u, C1) + + @Two_Product(acx, bcytail, s1, s0) + @Two_Product(acy, bcxtail, t1, t0) + @Two_Two_Diff(s1, s0, t1, t0, u3, u[3], u[2], u[1]) + u[4] = u3 + C2length = fast_expansion_sum_zeroelim(C1length, C1, 4, u, C2) + + @Two_Product(acxtail, bcytail, s1, s0) + @Two_Product(acytail, bcxtail, t1, t0) + @Two_Two_Diff(s1, s0, t1, t0, u3, u[3], u[2], u[1]) + u[4] = u3 + Dlength = fast_expansion_sum_zeroelim(C2length, C2, 4, u, D) + + return sign(D[Dlength]) +end + +function orientation(pa::Point{T}, pb::Point{T}, pc::Point{T}, + ccwerrboundA::Float64 = Devices.ccwerrboundA) where {T <: Float64Like} + detleft = (pa[1] - pc[1]) * (pb[2] - pc[2]) + detright = (pa[2] - pc[2]) * (pb[1] - pc[1]) + det = detleft - detright + + if detleft > zero(detleft) + if detright <= zero(detright) + return sign(det) + else + detsum = detleft + detright + end + elseif detleft < zero(detleft) + if detright >= zero(detright) + return sign(det) + else + detsum = -detleft - detright + end + else + return sign(det) + end + + errbound = ccwerrboundA * detsum + if (det >= errbound) || (-det >= errbound) + return sign(det) + end + + return orient2dadapt(pa, pb, pc, detsum) +end + +# h cannot be e or f. +function fast_expansion_sum_zeroelim(elen::Int, e, flen::Int, f, h) + enow = e[1] + fnow = f[1] + eindex = findex = 1 + if (fnow > enow) == (fnow > -enow) + Q = enow + eindex += 1 + enow = e[eindex] + else + Q = fnow + findex += 1 + fnow = f[findex] + end + hindex = 0 + if (eindex <= elen) && (findex <= flen) + if (fnow > enow) == (fnow > -enow) + @Fast_Two_Sum(enow, Q, Qnew, hh) + eindex += 1 + if eindex <= elen + enow = e[eindex] + end + else + @Fast_Two_Sum(fnow, Q, Qnew, hh) + findex += 1 + if findex <= flen + fnow = f[findex] + end + end + Q = Qnew + if hh != zero(hh) + hindex += 1 + h[hindex] = hh + end + while (eindex <= elen) && (findex <= flen) + if (fnow > enow) == (fnow > -enow) + @Two_Sum(Q, enow, Qnew, hh) + eindex += 1 + if eindex <= elen + enow = e[eindex] + end + else + @Two_Sum(Q, fnow, Qnew, hh) + findex += 1 + if findex <= flen + fnow = f[findex] + end + end + Q = Qnew + if hh != zero(hh) + hindex += 1 + h[hindex] = hh + end + end + end + while eindex <= elen + @Two_Sum(Q, enow, Qnew, hh) + eindex += 1 + if eindex <= elen + enow = e[eindex] + end + Q = Qnew + if hh != zero(hh) + hindex += 1 + h[hindex] = hh + end + end + while findex <= flen + @Two_Sum(Q, fnow, Qnew, hh) + findex += 1 + if findex <= flen + fnow = f[findex] + end + Q = Qnew + if hh != zero(hh) + hindex += 1 + h[hindex] = hh + end + end + if (Q != zero(Q)) || (hindex == 0) + hindex += 1 + h[hindex] = Q + end + return hindex +end diff --git a/src/render/compound.jl b/src/render/compound.jl index afaa837..09cc51a 100644 --- a/src/render/compound.jl +++ b/src/render/compound.jl @@ -9,10 +9,16 @@ function render!(c::Cell, f, len, s::Paths.CompoundStyle, meta::Meta; kwargs...) end end -function render!(c::Cell, seg::Paths.CompoundSegment, s::Paths.CompoundStyle, meta::Meta; +function render!(c::Cell, seg::Paths.CompoundSegment, sty::Paths.CompoundStyle, meta::Meta; kwargs...) - @assert length(seg.segments) == length(s.styles) - for (se,st) in zip(seg.segments, s.styles) - render!(c, se, st, meta; kwargs...) + if seg.tag == sty.tag + # Tagging mechanism: if some nodes have been simplified, but neither the segment + # nor the style have been swapped out, we'll just render as if we never simplified. + for (se, st) in zip(seg.segments, sty.styles) + render!(c, se, st, meta; kwargs...) + end + else + # Something has been done post-simplification, so we fall back to generic rendering + render!(c, seg, pathlength(seg), sty, meta; kwargs...) end end diff --git a/src/render/corners.jl b/src/render/corners.jl index f8d737a..988ab21 100644 --- a/src/render/corners.jl +++ b/src/render/corners.jl @@ -1,12 +1,16 @@ function render!(c::Cell, seg::Paths.Corner, ::Paths.SimpleTraceCorner, meta::Meta) sgn = ifelse(seg.α >= 0.0°, 1, -1) + + ext = seg.extent*tan(sgn*seg.α/2) + p0 = seg.p0 - ext*Point(cos(seg.α0), sin(seg.α0)) + ∠A = seg.α0+sgn*π/2 p = Point(cos(∠A),sin(∠A)) - p1 = seg.extent*p + seg.p0 - p2 = -seg.extent*p + seg.p0 - ex = 2*seg.extent*tan(sgn*seg.α/2) - p3 = p2 + ex*Point(cos(seg.α0),sin(seg.α0)) - p4 = p3 + ex*Point(cos(seg.α0+seg.α), sin(seg.α0+seg.α)) + + p1 = seg.extent*p + p0 + p2 = -seg.extent*p + p0 + p3 = p2 + 2ext*Point(cos(seg.α0),sin(seg.α0)) + p4 = p3 + 2ext*Point(cos(seg.α0+seg.α), sin(seg.α0+seg.α)) render!(c, Polygon([p1,p2,p3,p4]), Polygons.Plain(), meta) end diff --git a/src/render/paths.jl b/src/render/paths.jl index 19c2652..c44afd0 100644 --- a/src/render/paths.jl +++ b/src/render/paths.jl @@ -3,55 +3,16 @@ Render a path `p` to a cell `c`. """ function render!(c::Cell, p::Path{T}, meta::Meta=GDSMeta(); kwargs...) where {T} - inds = findall(x->isa(segment(x), Paths.Corner), nodes(p)) - segs = [] - - # Adjust the path so corners, when rendered with finite extent, - # are properly positioned. - # TODO: Add error checking for styles. - - for i in inds - cornernode = p[i] - prevseg = segment(previous(cornernode)) - nextseg = segment(next(cornernode)) - segs = [segs; prevseg; nextseg] - cornertweaks!(cornernode, prevseg, previous) - cornertweaks!(cornernode, nextseg, next) - end - !isempty(inds) && reconcile!(p) taper_inds = handle_generic_tapers!(p) - for node in p render!(c, segment(node), style(node), meta; kwargs...) end - - # Restore corner positions - for i in reverse(inds) - setsegment!(next(p[i]), pop!(segs)) - setsegment!(previous(p[i]), pop!(segs)) - end - !isempty(inds) && reconcile!(p) - restore_generic_tapers!(p, taper_inds) return c end -function cornertweaks!(cornernode, seg::Paths.Straight, which) - seg′ = copy(seg) - setsegment!(which(cornernode), seg′) - - α = segment(cornernode).α - ex = segment(cornernode).extent - sgn = ifelse(α >= 0.0°, 1, -1) - seg′.l -= ex*tan(sgn*α/2) -end - -cornertweak!(cornernode, seg::Paths.Segment) = - warn("corner was not sandwiched by straight segments. ", - "Rendering errors will result.") - function handle_generic_tapers!(p) # Adjust the path so generic tapers render correctly generic_taper_inds = findall(x->isa(style(x), Paths.Taper), nodes(p)) diff --git a/src/utils.jl b/src/utils.jl index c875fc1..205bf4b 100644 --- a/src/utils.jl +++ b/src/utils.jl @@ -39,7 +39,7 @@ adapted_grid(f, anchors::AbstractVector{T}; function assemble_grids(f, anchors::AbstractVector, max_recursions, max_change, rand_factor, grid_step) - dimension(f(anchors[1])) != dimension(max_change) && + dimension(eltype(f(anchors[1]))) != dimension(max_change) && throw(ArgumentError("max_change must have dimensions of f($(anchors[1]))")) g = i->begin # Want a range that specially excludes the last point; use linspace with npts. diff --git a/test/runtests.jl b/test/runtests.jl index b0a580f..0199ac0 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -493,6 +493,16 @@ end # end end +@testset "Intersections" begin + @test promote(Polygons.Line(p(0,0), p(1,1)), Polygons.Line(p(1.0,0.0), p(2.0,0.0))) === + (Polygons.Line(p(0.0,0.0), p(1.0,1.0)), Polygons.Line(p(1.0,0.0), p(2.0,0.0))) + @test promote( + Polygons.Line(p(0,0)μm2μm, p(1,1)μm2μm), + Polygons.Line(p(1.0,0.0)nm2μm, p(2.0,0.0)nm2μm)) === + (Polygons.Line(p(0.0,0.0)μm2μm, p(1.0,1.0)μm2μm), + Polygons.Line(p(0.001,0.0)μm2μm, p(0.002,0.0)μm2μm)) +end + @testset "Cell methods" begin # Setup nested cell refs c = Cell("main") @@ -729,7 +739,7 @@ end @test grid == [0.0μm, 25μm, 50μm, 75μm, 100μm] end -@testset "Style rendering" begin +@testset "Styles" begin @testset "NoRender" begin c = Cell("main") pa = Path(α0=24.31°) @@ -766,10 +776,10 @@ end p(-5.0nm, 5.0nm) ] @test points(c.elements[2]) == Point{typeof(1.0nm)}[ - p(34142.13562373095nm, 5850.793308457187nm), - p(34149.206691542815nm, 5857.864376269053nm), - p(34142.13562373095nm, 5864.9354440809175nm), - p(34135.06455591909nm, 5857.864376269053nm) + p(34142.13562373095nm, 5850.793308457185nm), + p(34149.206691542815nm, 5857.864376269051nm), + p(34142.13562373095nm, 5864.935444080916nm), + p(34135.06455591909nm, 5857.864376269051nm) ] @test points(c.elements[3]) == Point{typeof(1.0nm)}[ p(40005.0nm, 39995.0nm), @@ -798,10 +808,10 @@ end p(-4.99999999999997nm, 10505.0nm) ] @test points(c.elements[2]) == Point{typeof(1.0nm)}[ - p(26717.5144212722nm, 13275.414510915934nm), - p(26724.585489084064nm, 13282.485578727801nm), - p(26717.5144212722nm, 13289.556646539668nm), - p(26710.443353460334nm, 13282.485578727801nm) + p(26717.5144212722nm, 13275.414510915933nm), + p(26724.585489084064nm, 13282.4855787278nm), + p(26717.5144212722nm, 13289.556646539666nm), + p(26710.443353460334nm, 13282.4855787278nm) ] @test points(c.elements[3]) == Point{typeof(1.0nm)}[ p(29505.0nm, 39995.0nm), @@ -957,28 +967,28 @@ end turn!(pa, π/2, 50.0, Paths.CPW(10.0,6.0)) render!(c, pa, grid_step=500.0, max_change=90°) @test points(c.elements[1]) == Point{Float64}[ - p(6.79678973526781e-15, 11.0), - p(14.147528631528852, 13.656535200670929), - p(27.587463085789658, 22.433138000668393), - p(36.518103010286126, 36.30955981240446), + p(6.735557395310443e-16, 11.0), + p(14.14752863152884, 13.656535200670923), + p(27.58746308578965, 22.433138000668386), + p(36.518103010286126, 36.309559812404444), p(39.0, 50.0), p(45.0, 50.0), - p(42.13627270417631, 34.203338245082065), - p(31.83168817591114, 18.19208230846353), - p(16.324071497917906, 8.06523292385107), - p(6.429395695523605e-15, 5.0) + p(42.13627270417631, 34.20333824508206), + p(31.831688175911133, 18.192082308463522), + p(16.324071497917892, 8.065232923851063), + p(3.061616997868383e-16, 5.0) ] @test points(c.elements[2]) == Point{Float64}[ - p(6.429395695523605e-15, -5.0), - p(19.95164294189966, -1.2536042041820252), - p(38.90539665944695, 11.123656154788758), - p(51.49988886065992, 30.69296896621141), - p(55.0, 50.0), - p(61.0, 50.0), - p(57.1180585545501, 28.58674739888902), - p(43.14962174956843, 6.882600462583895), - p(22.128185808288713, -6.844906481001884), - p(6.79678973526781e-15, -11.0) + p(3.061616997868383e-16, -5.0), + p(19.951642941899642, -1.2536042041820332), + p(38.90539665944694, 11.123656154788751), + p(51.49988886065992, 30.692968966211403), + p(55.0, 50.0), + p(61.0, 50.0), + p(57.1180585545501, 28.586747398889013), + p(43.149621749568425, 6.882600462583888), + p(22.128185808288695, -6.8449064810018925), + p(6.735557395310443e-16, -11.0) ] c = Cell("main", nm) @@ -986,30 +996,49 @@ end turn!(pa, π/2, 50.0μm, Paths.CPW(10.0μm, 6.0μm)) render!(c, pa, grid_step=500.0μm, max_change=90°) @test points(c.elements[1]) == Point{typeof(1.0nm)}[ - p(6.79678973526781e-12nm, 11000.0nm), - p(14147.528631528852nm, 13656.535200670929nm), - p(27587.46308578966nm, 22433.138000668394nm), - p(36518.103010286126nm, 36309.559812404455nm), - p(39000.0nm, 50000.0nm), - p(45000.0nm, 50000.0nm), - p(42136.27270417631nm, 34203.338245082065nm), - p(31831.68817591114nm, 18192.08230846353nm), - p(16324.071497917907nm, 8065.23292385107nm), - p(6.429395695523605e-12nm, 5000.0nm) + (6.735557395310443e-13nm, 11000.0nm), + (14147.52863152884nm, 13656.535200670924nm), + (27587.46308578965nm, 22433.138000668387nm), + (36518.103010286126nm, 36309.55981240444nm), + (39000.0nm, 50000.0nm), + (45000.0nm, 50000.0nm), + (42136.27270417631nm, 34203.33824508206nm), + (31831.68817591113nm, 18192.082308463523nm), + (16324.071497917892nm, 8065.232923851063nm), + (3.061616997868383e-13nm, 5000.0nm) ] @test points(c.elements[2]) == Point{typeof(1.0nm)}[ - p(6.429395695523605e-12nm, -5000.0nm), - p(19951.64294189966nm, -1253.6042041820251nm), - p(38905.39665944695nm, 11123.656154788758nm), - p(51499.88886065992nm, 30692.96896621141nm), - p(55000.0nm, 50000.0nm), - p(61000.0nm, 50000.0nm), - p(57118.0585545501nm, 28586.74739888902nm), - p(43149.62174956843nm, 6882.600462583895nm), - p(22128.185808288712nm, -6844.906481001884nm), - p(6.79678973526781e-12nm, -11000.0nm) + p(3.061616997868383e-13nm, -5000.0nm), + p(19951.64294189964nm, -1253.6042041820333nm), + p(38905.396659446946nm, 11123.65615478875nm), + p(51499.88886065992nm, 30692.968966211403nm), + p(55000.0nm, 50000.0nm), + p(61000.0nm, 50000.0nm), + p(57118.0585545501nm, 28586.747398889012nm), + p(43149.621749568425nm, 6882.600462583888nm), + p(22128.185808288694nm, -6844.906481001893nm), + p(6.735557395310443e-13nm, -11000.0nm) ] + pa = Path(μm) + turn!(pa, π/2, 50.0μm, Paths.CPW(10.0μm, 6.0μm)) + + pa2 = split(pa[1], 50.0μm * 30°) + let s1 = style(pa2[1]), s2 = style(pa2[2]) + @test Paths.trace(s1, 0μm) == 10.0μm + @test Paths.trace(s1, 50.0μm * 30°) == 10.0μm + @test Paths.trace(s2, 0μm) == 10.0μm + @test Paths.trace(s2, 50.0μm * 60°) == 10.0μm + @test Paths.gap(s1, 0μm) == 6.0μm + @test Paths.gap(s1, 50.0μm * 30°) == 6.0μm + @test Paths.gap(s2, 0μm) == 6.0μm + @test Paths.gap(s2, 50.0μm * 60°) == 6.0μm + end + let s1 = segment(pa2[1]), s2 = segment(pa2[2]) + @test p0(s1) == Point(0,0)μm + @test p1(s1) == p0(s2) ≈ Point(50.0*sin(30°), 50*(1-cos(30°)))μm + @test p1(s2) ≈ Point(50,50)μm + end end @testset "Straight, TaperTrace" begin @@ -1018,11 +1047,29 @@ end straight!(pa, 50.0μm, Paths.TaperTrace(10.0μm, 6.0μm)) render!(c, pa) @test points(c.elements[1]) ≈ Point{typeof(1.0nm)}[ - p(0.0nm, 5000.0nm), - p(50000.0nm, 3000.0nm), + p(0.0nm, 5000.0nm), + p(50000.0nm, 3000.0nm), p(50000.0nm, -3000.0nm), - p(0.0nm, -5000.0nm) + p(0.0nm, -5000.0nm) ] + + # length not yet specified + @test_throws AssertionError split(Paths.TaperTrace(10.0μm, 6.0μm), 10μm) + + pa2 = split(pa[1], 10μm) + let s1 = style(pa2[1]), s2 = style(pa2[2]) + @test Paths.width(s1, 0μm) ≈ 10.0μm + @test Paths.width(s1, 10μm) ≈ 9.2μm + @test s1.length == 10μm + @test Paths.width(s2, 0μm) ≈ 9.2μm + @test Paths.width(s2, 40μm) ≈ 6.0μm + @test s2.length == 40μm + end + let s1 = segment(pa2[1]), s2 = segment(pa2[2]) + @test p0(s1) == Point(0,0)μm + @test p1(s1) == p0(s2) == Point(10,0)μm + @test p1(s2) == Point(50,0)μm + end end @testset "Straight, TaperCPW" begin @@ -1031,17 +1078,33 @@ end straight!(pa, 50.0μm, Paths.TaperCPW(10.0μm, 6.0μm, 8.0μm, 2.0μm)) render!(c, pa) @test points(c.elements[1]) ≈ Point{typeof(1.0nm)}[ - p(0.0nm, 11000.0nm), - p(50000.0nm, 6000.0nm), - p(50000.0nm, 4000.0nm), - p(0.0nm, 5000.0nm) + p(0.0nm, 11000.0nm), + p(50000.0nm, 6000.0nm), + p(50000.0nm, 4000.0nm), + p(0.0nm, 5000.0nm) ] @test points(c.elements[2]) ≈ Point{typeof(1.0nm)}[ - p(0.0nm, -5000.0nm), - p(50000.0nm, -4000.0nm), - p(50000.0nm, -6000.0nm), - p(0.0nm, -11000.0nm) + p(0.0nm, -5000.0nm), + p(50000.0nm, -4000.0nm), + p(50000.0nm, -6000.0nm), + p(0.0nm, -11000.0nm) ] + + @test_throws AssertionError split(Paths.TaperCPW(10.0μm, 6.0μm, 8.0μm, 2.0μm), 10μm) + + pa2 = split(pa[1], 10μm) + let s1 = style(pa2[1]), s2 = style(pa2[2]) + @test Paths.trace(s1, 0μm) ≈ 10.0μm + @test Paths.trace(s1, 10μm) ≈ 9.6μm + @test Paths.gap(s1, 0μm) ≈ 6.0μm + @test Paths.gap(s1, 10μm) ≈ 5.2μm + @test s1.length == 10μm + @test Paths.trace(s2, 0μm) ≈ 9.6μm + @test Paths.trace(s2, 40μm) ≈ 8.0μm + @test Paths.gap(s2, 0μm) ≈ 5.2μm + @test Paths.gap(s2, 40μm) ≈ 2.0μm + @test s2.length == 40μm + end end @testset "CompoundSegment" begin @@ -1082,6 +1145,21 @@ end straight!(pa, 20μm, Paths.Trace(15μm)) straight!(pa, 20μm, Paths.Trace(20μm)) simplify!(pa) + + pa2 = split(pa[1], 20μm) + @test length(pa2) == 2 + @test length(segment(pa2[1]).segments) == 1 + @test p1(segment(pa2[1])) == p0(segment(pa2[2])) == Point(20,0)μm + @test p1(segment(pa2[2])) == Point(60,0)μm + @test length(segment(pa2[2]).segments) == 2 + + pa2 = split(pa[1], 30μm) + @test length(pa2) == 2 + @test length(segment(pa2[1]).segments) == 2 + @test p1(segment(pa2[1])) == p0(segment(pa2[2])) == Point(30,0)μm + @test p1(segment(pa2[2])) == Point(60,0)μm + @test length(segment(pa2[2]).segments) == 2 + setsegment!(pa[1], Paths.Straight(120.0μm, p(0.0μm, 0.0μm), 0.0)) render!(c, pa) @test lowerleft(bounds(c.elements[1])) ≈ Point(0μm, -5μm)