From e98fa1fddedf21b609474445c3c421af82013f90 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Sat, 17 Oct 2026 20:01:15 +0200 Subject: [PATCH 01/20] Bump MemoryViews --- Project.toml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index d64435c..76ccc0f 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "PairwiseMappingFormat" uuid = "30ee2946-b067-41c0-bed2-0900b78dbe59" -version = "0.1.2" +version = "0.1.3" authors = ["Jakob Nybo Nissen "] [deps] @@ -12,7 +12,7 @@ XAMAuxData = "e99d641e-1821-45d7-9150-ecb7bf333fe1" [compat] Aqua = "0.8.7" FormatSpecimens = "1.3" -MemoryViews = "0.2, 0.3, 0.4" +MemoryViews = "0.2, 0.3, 0.4, 0.5" PrecompileTools = "1.2.1" StringViews = "1.3.3" Test = "1.11" From ad5418518abc29fbcc2c30b5b6b216d6ef6d5aa9 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Sat, 17 Oct 2026 20:30:45 +0200 Subject: [PATCH 02/20] WIP: v0.2.0 --- CHANGELOG.md | 14 ++++++++++++++ Project.toml | 10 ++++++---- src/PairwiseMappingFormat.jl | 3 ++- test/runtests.jl | 7 +++++-- 4 files changed, 27 insertions(+), 7 deletions(-) create mode 100644 CHANGELOG.md diff --git a/CHANGELOG.md b/CHANGELOG.md new file mode 100644 index 0000000..c4e8247 --- /dev/null +++ b/CHANGELOG.md @@ -0,0 +1,14 @@ +# Changelog +All notable changes to this project will be documented in this file. + +The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.1.0/), +and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html). + +## v0.2.0 +### Changed +- **Breaking:** Allow XAMAuxData 0.2. Since `aux_data` returns a `XAMAuxData.SAM.Auxiliary`, + the breaking changes of XAMAuxData 0.2 apply to auxiliary data of PAF records. + For example, `H` values with lowercase hex digits are now loaded as an `Error`, + and `length` and `haskey` throw on malformed auxiliary data. + See the migration guide in the documentation of XAMAuxData.jl. +- Allow MemoryViews 0.5. diff --git a/Project.toml b/Project.toml index 76ccc0f..2b50890 100644 --- a/Project.toml +++ b/Project.toml @@ -1,9 +1,10 @@ name = "PairwiseMappingFormat" uuid = "30ee2946-b067-41c0-bed2-0900b78dbe59" -version = "0.1.3" -authors = ["Jakob Nybo Nissen "] +version = "0.2.0" +authors = ["Jakob Nybo Andersen "] [deps] +BufferIO = "ff366017-094f-4144-bd70-a828a0d67152" MemoryViews = "a791c907-b98b-4e44-8f4d-e4c2362c6b2f" PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a" StringViews = "354b36f9-a18e-4713-926e-db85100087ba" @@ -11,12 +12,13 @@ XAMAuxData = "e99d641e-1821-45d7-9150-ecb7bf333fe1" [compat] Aqua = "0.8.7" +BufferIO = "0.2.7" FormatSpecimens = "1.3" MemoryViews = "0.2, 0.3, 0.4, 0.5" PrecompileTools = "1.2.1" -StringViews = "1.3.3" +StringViews = "2.0.0" Test = "1.11" -XAMAuxData = "0.1" +XAMAuxData = "0.2" julia = "1.11" [extras] diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index b2980d5..2a57f49 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -25,6 +25,7 @@ end using MemoryViews: MemoryView, ImmutableMemoryView using XAMAuxData.SAM: Auxiliary +using BufferIO: ExpandReader """ PAFRecord(buffer_size::Int) @@ -62,7 +63,7 @@ The following properties may be used: * `matches::Int` gives the number of nucleotides / residues which match (i.e. are equal) in the alignment. * `alnlen::Int` gives the length of the alignment -* `mapq::Union{Int, Nothing}` gives the mapping quality, or `nothing` if +* `mapq::Union{Int, Nothing}` gives the mapping quality, or `nothing` if this information is unavailable. * `is_rc::Union{Bool,Nothing}` tells if the query and target strands are reverse-complement relative to each other. Is `nothing` if the record is unmapped. diff --git a/test/runtests.jl b/test/runtests.jl index 6d1b5bb..51cf016 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -7,6 +7,7 @@ using FormatSpecimens using MemoryViews using Aqua +using XAMAuxData: XAMAuxData using XAMAuxData.SAM: AuxTag function cmp_aux(record::PAFRecord, other) @@ -15,7 +16,7 @@ function cmp_aux(record::PAFRecord, other) @test length(aux) == length(other) for ((kaux, vaux), (ko, vo)) in zip(aux, other) @test kaux == AuxTag(ko) - if vaux isa Union{Real, AbstractString, Char} + if vaux isa Union{Real, AbstractString, Char, XAMAuxData.Error} @test vaux == vo elseif vaux isa AbstractVector @test eltype(vaux) == eltype(vo) @@ -91,7 +92,9 @@ end end @test rec1.qname == "中国科学院大学" @test rec1.qlen == typemax(Int32) - cmp_aux(rec1, ["tp" => 'S', "cm" => 29990, "s1" => 301142, "dv" => -13.2442f-18, "rl" => 0, "xX" => [0xea, 0x19, 0xfe, 0x1b], "K1" => Int8[15, 9, 22, 127, -12, -128, 0]]) + # XAMAuxData 0.2 rejects lowercase hex digits + xX = pkgversion(XAMAuxData) < v"0.2" ? [0xea, 0x19, 0xfe, 0x1b] : XAMAuxData.Errors.InvalidHex + cmp_aux(rec1, ["tp" => 'S', "cm" => 29990, "s1" => 301142, "dv" => -13.2442f-18, "rl" => 0, "xX" => xX, "K1" => Int8[15, 9, 22, 127, -12, -128, 0]]) cmp_aux(rec2, ["dv" => 1.243f-4, "rl" => 38, "kk" => "!=[]! And then~~ "]) cmp_aux(rec3, ["dv" => +17.0f21, "rl" => 38]) From ef41b908b72901d6b7bbd13cc4935e18333ae1e4 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Sat, 17 Oct 2026 21:21:31 +0200 Subject: [PATCH 03/20] Bump StringViews, XAMAuxData --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 2b50890..265d222 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "PairwiseMappingFormat" uuid = "30ee2946-b067-41c0-bed2-0900b78dbe59" -version = "0.2.0" +version = "0.1.2" authors = ["Jakob Nybo Andersen "] [deps] From f843cfc48df41320ee5cbd269ffe52f4cdf56b35 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Sat, 17 Oct 2026 21:10:14 +0200 Subject: [PATCH 04/20] Back reader with AbstractBufReader --- CHANGELOG.md | 18 +++ Project.toml | 3 +- docs/src/index.md | 17 +-- src/PairwiseMappingFormat.jl | 257 ++++++++++++----------------------- test/runtests.jl | 27 ++++ 5 files changed, 140 insertions(+), 182 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index c4e8247..b932232 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -12,3 +12,21 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 and `length` and `haskey` throw on malformed auxiliary data. See the migration guide in the documentation of XAMAuxData.jl. - Allow MemoryViews 0.5. +- **Breaking:** `PAFReader` now reads from a `BufferIO.AbstractBufReader`, and its type + parameter is the type of the buffered reader. An `IO` passed to `PAFReader` is + wrapped in a `BufferIO.BufReader`, and so becomes `PAFReader{BufReader{T}}`. + Reading from an `AbstractBufReader` avoids copying the data into the reader's own buffer. +- **Breaking:** A `PAFRecord` may now be a view into the memory it was parsed from, instead + of owning a copy of the data: + - With `copy=false`, the record iterated by `PAFReader` and returned by `try_next!` + is a view into the reader's buffer, and is invalidated when the reader is advanced. + - `try_parse` and `parse` return a view into the input, so mutating the input + also mutates the record. + - `copy(record)` returns a record that owns its data. +- **Breaking:** `aux_data` now returns an immutable `Auxiliary{ImmutableMemoryView{UInt8}}`, + so the auxiliary data of records can no longer be modified. +- **Breaking:** Reading from a closed `PAFReader` now throws, instead of + behaving as if the reader is empty. + +### Removed +- **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. diff --git a/Project.toml b/Project.toml index 265d222..d607b46 100644 --- a/Project.toml +++ b/Project.toml @@ -23,6 +23,7 @@ julia = "1.11" [extras] Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" +BufferIO = "ff366017-094f-4144-bd70-a828a0d67152" FormatSpecimens = "3372ea36-2a1a-11e9-3eb7-996970b6ffbd" MemoryViews = "a791c907-b98b-4e44-8f4d-e4c2362c6b2f" StringViews = "354b36f9-a18e-4713-926e-db85100087ba" @@ -30,4 +31,4 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" XAMAuxData = "e99d641e-1821-45d7-9150-ecb7bf333fe1" [targets] -test = ["Test", "FormatSpecimens", "MemoryViews", "StringViews", "XAMAuxData", "Aqua"] +test = ["Test", "FormatSpecimens", "MemoryViews", "StringViews", "XAMAuxData", "Aqua", "BufferIO"] diff --git a/docs/src/index.md b/docs/src/index.md index 695211d..8a7f6a0 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -16,14 +16,15 @@ which is a simple tab-delimited format used by e.g. minimap2 and strobealign to store records representing pairwise alignment between biological sequences. ## Reader -The [`PAFReader`](@ref) type wraps an `IO`, and is an iterator of [`PAFRecord`](@ref) objects: +The [`PAFReader`](@ref) type wraps an `IO` or a `BufferIO.AbstractBufReader`, +and is an iterator of [`PAFRecord`](@ref) objects: ```jldoctest; output = false reader = PAFReader(open(path_to_paf)) records = collect(reader) +@assert isempty(reader) close(reader) -@assert isempty(reader) @assert typeof(records) == Vector{PAFRecord} # output @@ -47,9 +48,10 @@ julia> PAFReader(open(path_to_paf)) do reader The [`PAFReader`](@ref) constructor takes an optional keyword `copy`, which defaults to `true`. If it is `false`, the record will iterate the _same_ [`PAFRecord`](@ref) object, overwriting -it in-place. -This reduces allocations and gives a slight speed increase, but may result in bugs -if records, or references of records of previous iterations are stored and +it in-place. This record is a view into the reader's buffer, which is only valid until the +reader is advanced. +This reduces allocations and copying and gives a slight speed increase, but may result +in bugs if records, or references of records of previous iterations are stored and unexpectedly overwritten: ```jldoctest @@ -266,11 +268,6 @@ julia> aux["km"] julia> haskey(aux, "cm") true - -julia> aux["AB"] = 5.5; - -julia> haskey(aux, "AB") -true ``` ## Misc info diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 2a57f49..4bbfd96 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -22,17 +22,21 @@ using PrecompileTools: @compile_workload, @setup_workload, @recompile_invalidati using StringViews: StringView end -using MemoryViews: MemoryView, ImmutableMemoryView +using MemoryViews: ImmutableMemoryView using XAMAuxData.SAM: Auxiliary -using BufferIO: ExpandReader +using BufferIO: AbstractBufReader, BufReader, ExpandReader, get_buffer, fill_buffer, consume """ - PAFRecord(buffer_size::Int) + PAFRecord Mutable type representing a single PAF line. The content of the record is accessed through its properties. +A `PAFRecord` does not necessarily own its data, but may be a view into the +memory it was parsed from, e.g. the buffer of a [`PAFReader`](@ref) with `copy=false`. +Use `copy` to obtain a record that owns its data. + # Examples ```jldoctest julia> record = PAFReader(first, open(path_to_paf)); @@ -69,11 +73,13 @@ The following properties may be used: relative to each other. Is `nothing` if the record is unmapped. """ mutable struct PAFRecord - # Data contains query name, subject name, AUX data, - # without delimiters and with no trailing whitespace. - const data::Vector{UInt8} - qname_len::Int32 - tname_len::Int32 + # The line the record was parsed from, without the trailing newline. + # May be a view into memory not owned by the record, e.g. a reader's buffer. + data::ImmutableMemoryView{UInt8} + qname_len::Int32 # query name begins at index 1 + tname_start::Int32 + tname_len::Int32 # always zero when unmapped + aux_start::Int32 # lastindex(data) + 1 if there is no aux data qlen::Int32 # always >= 1 # Note: qstart/qend (and tstart/tend) are zero-indexed semi-open interval # in PAF format, but stored in this struct one-based. @@ -150,9 +156,9 @@ function Base.propertynames(::PAFRecord) ) end -function PAFRecord(size::Int = 0) - data = Vector{UInt8}(undef, max(size, 0)) - return PAFRecord(data, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0xff, 0x00) +function PAFRecord() + data = ImmutableMemoryView(Memory{UInt8}()) + return PAFRecord(data, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0xff, 0x00) end function Base.show(io::IO, ::MIME"text/plain", record::PAFRecord) @@ -188,17 +194,15 @@ end # Return StringView to not allocate function qname(record::PAFRecord) - return StringView( - ImmutableMemoryView(getfield(record, :data))[1:(getfield(record, :qname_len))], - ) + return StringView(getfield(record, :data)[1:(getfield(record, :qname_len))]) end # Return StringView to not allocate function tname(record::PAFRecord)::Union{StringView, Nothing} return if is_mapped(record) - ql = getfield(record, :qname_len) - span = (ql + 1):(ql + getfield(record, :tname_len)) - StringView(ImmutableMemoryView(getfield(record, :data))[span]) + start = getfield(record, :tname_start) + span = start:(start + getfield(record, :tname_len) - 1) + StringView(getfield(record, :data)[span]) else nothing end @@ -281,17 +285,12 @@ julia> length(aux) julia> aux["cm"] 29990 -julia> aux["k1"] = "add a new aux string like this"; - julia> haskey(aux, "k1") -true +false ``` """ function aux_data(record::PAFRecord) - return Auxiliary( - getfield(record, :data), - getfield(record, :qname_len) + getfield(record, :tname_len) + 1, - ) + return Auxiliary(getfield(record, :data)[getfield(record, :aux_start):end]) end """ @@ -410,6 +409,11 @@ end Try parsing `x` into a `PAFRecord`. `x` may be any type that implements `MemoryView`, such as `String` or `Vector{UInt8}`. +!!! warning + The resulting record is a view into the memory of `x`, and does not copy it. + If `x` is mutated, the record is mutated too, and may become invalid. + Use `copy` on the record to obtain a record that owns its data. + # Examples ```jldoctest julia> valid_line = open(readline, path_to_paf); @@ -424,9 +428,7 @@ PairwiseMappingFormat.ParserException See also: [`PAFRecord`](@ref), [`try_next!`](@ref) """ try_parse(x) = try_parse(ImmutableMemoryView(x)) -# At least 19 bytes of input are stored in the fixed fields, so -# we need at most 19 fewer bytes in the data field of the vector. -try_parse(x::ImmutableMemoryView{UInt8}) = parse_line!(PAFRecord(length(x) - 19), x) +try_parse(x::ImmutableMemoryView{UInt8}) = parse_line!(PAFRecord(), x) function Base.parse(::Type{PAFRecord}, s::Union{String, SubString{String}}) y = try_parse(s) @@ -437,8 +439,9 @@ function parse_line!( record::PAFRecord, mem::ImmutableMemoryView{UInt8}, )::Union{PAFRecord, ParserException} - if lastindex(mem) > typemax(Int32) - error("PairwiseMappingFormat.jl can't handle lines longer than 2147483647 bytes") + # Note: The aux data may begin at lastindex(mem) + 1, which must fit in an Int32 + if lastindex(mem) >= typemax(Int32) + error("PairwiseMappingFormat.jl can't handle lines longer than 2147483646 bytes") end # Read fields incrementally. `i` is the index to read from, beginning at 1. # The parse_x functions will return the next index to begin from (after the tab) @@ -484,24 +487,14 @@ function parse_line!( mapq % UInt8 end - # Copy data to the .data field - query name, target name, AUX data - auxlen = lastindex(mem) - i + 1 - filled = length(qname) + length(tname) + auxlen - data = getfield(record, :data) - length(data) == filled || resize!(data, filled) - doff = 1 - dataview = MemoryView(data) - unsafe_copyto!(dataview, mem[qname]) - doff += length(qname) - unsafe_copyto!(dataview[doff:end], mem[tname]) - doff += length(tname) - iszero(auxlen) || unsafe_copyto!(dataview[doff:end], mem[i:end]) - - # Fill in fields of the struct + # Fill in fields of the struct. The aux data begins at `i`. # Note: PAF uses zero-based semi-open indexing like e.g. Python, # so compensate for that in qstart and qend. + record.data = mem record.qname_len = length(qname) % Int32 + record.tname_start = first(tname) % Int32 record.tname_len = length(tname) % Int32 + record.aux_start = i % Int32 record.qlen = qlen record.qstart = qstart + one(qstart) record.qend = qend @@ -536,21 +529,14 @@ function finish_unmapped!( end end i = findnext(==(UInt8('\t')), mem, i) - aux_end = lastindex(mem) - aux_start = isnothing(i) ? aux_end + 1 : i + 1 - aux = aux_start:aux_end - - # Copy data - filled = length(qname) + length(aux) - data = getfield(record, :data) - dataview = MemoryView(data) - length(data) == filled || resize!(data, filled) - unsafe_copyto!(dataview, mem[qname]) - unsafe_copyto!(dataview[(length(qname) + 1):end], mem[aux]) + aux_start = isnothing(i) ? lastindex(mem) + 1 : i + 1 # Fill in fields + record.data = mem record.qname_len = length(qname) % Int32 + record.tname_start = 1 record.tname_len = 0 + record.aux_start = aux_start % Int32 record.qlen = qlen record.qstart = 0 record.qend = 0 @@ -615,21 +601,25 @@ function parse_int( end """ - PAFReader(io::IO; copy::Bool=true) + PAFReader(io::Union{IO, AbstractBufReader}; copy::Bool=true) Construct a `PAFReader`, an iterator of [`PAFRecord`](@ref) read from `io`. Iterating this object returns a `PAFRecord` for each valid line of input, and throws a [`ParserException`](@ref) when an invalid record is found. For efficiency, the auxiliary fields are not validated until they are accessed. -If `copy` is false, the *same* record will be returned on each iteration, -with its content being overwritten. This removes a few allocations per iteration, -but may cause problems if records of old iterations are stored. +If `io` is an `IO`, it is wrapped in a `BufferIO.BufReader` with a buffer size +given by the keyword `buf_size`. + +If `copy` is false, the *same* record will be returned on each iteration. +This record is a view into the reader's buffer, and is only valid until the reader is +advanced. This removes a few allocations and copies per iteration, but iterating the +reader will overwrite old records. # Examples ```jldoctest julia> reader = PAFReader(open(path_to_paf)); typeof(reader) -PAFReader{IOStream} +PAFReader{BufferIO.BufReader{IOStream}} julia> typeof(first(reader)) PAFRecord @@ -654,31 +644,31 @@ true See also: [`PAFRecord`](@ref), [`try_next!`](@ref) """ -mutable struct PAFReader{I <: IO} - const io::I +mutable struct PAFReader{I <: AbstractBufReader} + const io::ExpandReader{I} const record::PAFRecord - mem::Memory{UInt8} - # Index of first unused byte in buffer. Equal to filled + 1 - # if the buffer has no usable data. - index::Int32 - # mem[1:filled] contain usable data - filled::Int32 + # Number of bytes of the last parsed line, including its newline, which have not + # yet been consumed from `io`. Consuming is deferred to the next call to `try_next!`, + # such that the record, which is a view into the buffer of `io`, stays valid. + to_consume::Int line::Int # Whether to iterate copies of the record stored in this struct copy::Bool end +function PAFReader(io::AbstractBufReader; copy::Bool = true) + return PAFReader{typeof(io)}(ExpandReader(io), PAFRecord(), 0, 1, copy) +end + function PAFReader(io::IO; buf_size::Int = 2^16, copy::Bool = true) - mem = Memory{UInt8}(undef, max(buf_size, 16)) - rec = PAFRecord(512) - return PAFReader{typeof(io)}(io, rec, mem, 1, 0, 1, copy) + return PAFReader(BufReader(io, max(buf_size, 16)); copy) end Base.IteratorSize(::Type{<:PAFReader}) = Base.SizeUnknown() Base.eltype(::Type{<:PAFReader}) = PAFRecord Base.close(reader::PAFReader) = close(reader.io) # TODO: Docstring -function PAFReader(f::Function, io::IO; kwargs...) +function PAFReader(f::Function, io::Union{IO, AbstractBufReader}; kwargs...) reader = PAFReader(io; kwargs...) return try f(reader) @@ -691,7 +681,9 @@ function Base.copy(record::PAFRecord) return PAFRecord( copy(getfield(record, :data)), getfield(record, :qname_len), + getfield(record, :tname_start), getfield(record, :tname_len), + getfield(record, :aux_start), getfield(record, :qlen), getfield(record, :qstart), getfield(record, :qend), @@ -717,80 +709,6 @@ function Base.iterate(reader::PAFReader, (::Nothing) = nothing) end end -# This is never called during normal operations, and only if a single line -# does not fit in the default buffer, which is highly unlikely. -function reallocate!(reader::PAFReader) - mem = reader.mem - newlen = length(mem) * 2 - if newlen >= typemax(Int32) - error("PAFReader encountered a line with more than 2147483646 bytes.") - end - new = Memory{UInt8}(undef, newlen) - unsafe_copyto!(new, 1, mem, 1, reader.filled) - reader.mem = new - return reader -end - -# Shift the data inside the buffer such that the usable -# data begins at index 1, thus making room for more -# data at the end of the buffer -function shift_buffer!(reader::PAFReader) - isone(reader.index) && return reader - mem = reader.mem - len = reader.filled - reader.index + 1 - copyto!(mem, 1, mem, reader.index, len) - reader.index = 1 - reader.filled = len - return reader -end - -# Read data into the buffer. Returns zero if and only if underlying -# stream is EOF (or buggy). -function read_buffer!(reader::PAFReader)::Int - mem = reader.mem - remaining = length(mem) - reader.filled - # If no space, try to shift data in the buffer to make room - if iszero(remaining) - shift_buffer!(reader) - remaining = length(mem) - reader.filled - # If still no space, the buffer is truly full, and we need to - # resize it - if iszero(remaining) - reallocate!(reader) - mem = reader.mem # it was change in reallocate! - remaining = length(mem) - reader.filled - end - end - - @assert !iszero(remaining) - # This part annoys me a lot. Julia has no efficient API to read bytes - # into a byte buffer. How on Earth can that be true? - # * read! will throw EOF if the memory is too short, and there is no - # way to query how many bytes are available - # * readbytes! is not optimised for Memory - # * unsafe_read has almost no promised semantics so practically can't be used. - # The approach here uses unsafe_read in a manner similar to TranscodingStreams.jl. - # It uses bytesavailable to check how many bytes to read and avoid buffer overflows etc, - # and when there are zero bytes available, it reads a single byte, hoping that - # more bytes will be available. - n_available = bytesavailable(reader.io) - if iszero(n_available) && !eof(reader.io) - # Read one byte and hope more bytes become available - byte = read(reader.io, UInt8) - @inbounds mem[reader.filled + 1] = byte - reader.filled += 1 - remaining -= 1 - n_available = bytesavailable(reader.io) - end - # Read with unsafe_read - n_bytes = min(n_available, remaining) - if !iszero(n_bytes) - GC.@preserve mem unsafe_read(reader.io, pointer(mem, reader.filled + 1), n_bytes) - reader.filled += n_bytes - end - return n_bytes -end - """ try_next!(reader::PAFReader) -> Union{Nothing, PAFRecord, ParserException} @@ -803,6 +721,9 @@ Try to read a record from the [`PAFReader`](@ref) and advance the reader. This m Even if the reader itself is not advanced, it may still fill its internal buffer, advancing the `IO` object it wraps. +The returned record is the reader's own record, which is a view into the reader's buffer. +It is overwritten, and its data invalidated, by the next call to `try_next!`. + # Examples ```jldoctest julia> reader = PAFReader(open(path_to_paf)); @@ -826,49 +747,43 @@ true ``` """ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} - # n_scanned: The number of bytes in buffer we have already memchrd searching for a newline - n_scanned = 0 - newline = findfirst( - ==(0x0a), - ImmutableMemoryView(reader.mem)[(reader.index + n_scanned):(reader.filled % Int)], - ) - n_bytes_read = 0 + io = reader.io + # Consume the previous line, now that the previous record is invalidated. + if !iszero(reader.to_consume) + consume(io, reader.to_consume) + reader.to_consume = 0 + end + buffer = get_buffer(io) + newline = findfirst(==(0x0a), buffer) # Keep reading until we find a newline or reach EOF while newline === nothing - n_scanned += reader.filled - reader.index - n_scanned + 1 - n_bytes_read = read_buffer!(reader) + # The number of bytes in buffer we have already memchrd searching for a newline + n_scanned = length(buffer) # This signals EOF - we might still have data in the buffer, if the underlying stream # ends with a record that doesn't end in a newline - iszero(n_bytes_read) && break - newline = findfirst( - ==(0x0a), - ImmutableMemoryView(reader.mem)[(reader.index + n_scanned):(reader.filled % Int)], - ) + iszero(fill_buffer(io)) && break + buffer = get_buffer(io) + newline = findnext(==(0x0a), buffer, n_scanned + 1) end - # If a newline was found, the current value of newline was relative to where - # the last search began. So, update it to account for the bytes searched in - # previous iterations of search for tab - newline === nothing || (newline += reader.index + n_scanned - 1) - mem = reader.mem has_windows_newline = false # Last index before \n, and before \r if the newline is \r\n as on Windows last_index = if newline === nothing - reader.filled + length(buffer) else if isone(newline) 0 else - has_windows_newline = @inbounds(mem[newline - 1]) == UInt8('\r') + has_windows_newline = @inbounds(buffer[newline - 1]) == UInt8('\r') newline - (1 + has_windows_newline) end end - parse_mem = ImmutableMemoryView(mem)[reader.index:last_index] + parse_mem = buffer[1:last_index] # If there is no memory to parse (and we have already verified that # there are no more bytes in the IO), we end if isempty(parse_mem) line = reader.line col = 0 - @inbounds for b in ImmutableMemoryView(mem)[reader.index:reader.filled] + @inbounds for b in buffer col += 1 line += b == UInt8('\n') col *= b != UInt8('\n') @@ -884,9 +799,9 @@ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException return ParserException(res.kind, getfield(res, :index_in_line), reader.line) # Else, we update the reader. This means `next!` will continuously return # an error if an error is found. - # If the newline was found, we compensate for the newline by skipping that byte in - # the next iteration, hence the extra `newline !== nothing`. - reader.index = last_index + 1 + has_windows_newline + (newline !== nothing) + # If the newline was found, we compensate for the newline by also consuming that byte, + # hence the extra `newline !== nothing`. + reader.to_consume = last_index + has_windows_newline + (newline !== nothing) reader.line += 1 return res end diff --git a/test/runtests.jl b/test/runtests.jl index 51cf016..38da693 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -6,6 +6,7 @@ using Test using FormatSpecimens using MemoryViews using Aqua +using BufferIO: CursorReader using XAMAuxData: XAMAuxData using XAMAuxData.SAM: AuxTag @@ -196,6 +197,32 @@ end end end +@testset "AbstractBufReader" begin + good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" + other = replace(good, "my_qname" => "other") + data = join([good, other, good], '\n') + + # Records iterated with copy=true do not alias the reader's buffer + records = PAFReader(collect, CursorReader(data)) + @test [r.qname for r in records] == ["my_qname", "other", "my_qname"] + + # With copy=false, the record is valid until the reader is advanced + PAFReader(CursorReader(data); copy = false) do reader + rec = try_next!(reader) + @test rec.qname == "my_qname" + @test try_next!(reader) === rec + @test rec.qname == "other" + end + + # Parsing a vector aliases it, but copying the record does not + v = collect(codeunits(good)) + rec = try_parse(v) + copied = copy(rec) + v[1:2] .= UInt8('X') + @test rec.qname == "XX_qname" + @test copied.qname == "my_qname" +end + # We assume that aux data is correctly tested in its own package, here we just test # that the aux data has the right span and byte content. function test_aux_content(record_data::String, target_aux::String) From 9283bd1fb67f4f74c3897ee706222492bf11e5a6 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Sun, 18 Oct 2026 20:48:04 +0200 Subject: [PATCH 05/20] Fix bugs --- .gitignore | 2 ++ src/PairwiseMappingFormat.jl | 48 ++++++++++++++++++++++++------------ test/runtests.jl | 42 ++++++++++++++++++++++++++++++- 3 files changed, 75 insertions(+), 17 deletions(-) diff --git a/.gitignore b/.gitignore index db690aa..d9d8f9b 100644 --- a/.gitignore +++ b/.gitignore @@ -28,3 +28,5 @@ Manifest.toml TODO.md *.code-workspace + +/REVIEW.md diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 4bbfd96..32c3f5c 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -320,6 +320,7 @@ module Errors EmptyInteger InvalidZero PositionOutOfBounds + LineTooLong end end # module Errors @@ -389,6 +390,8 @@ function Base.showerror(io::IO, err::ParserException) "Parsed an integer to zero which cannot meaningfully be zero" elseif kind == Errors.PositionOutOfBounds "Query or target start or end position is longer than query/target length" + elseif kind == Errors.LineTooLong + "Line is too long. PairwiseMappingFormat.jl can't handle lines of 2147483647 bytes or longer" end print(buf, s) return print(io, String(take!(buf))) @@ -440,15 +443,19 @@ function parse_line!( mem::ImmutableMemoryView{UInt8}, )::Union{PAFRecord, ParserException} # Note: The aux data may begin at lastindex(mem) + 1, which must fit in an Int32 - if lastindex(mem) >= typemax(Int32) - error("PairwiseMappingFormat.jl can't handle lines longer than 2147483646 bytes") - end + lastindex(mem) >= typemax(Int32) && + return ParserException(typemax(Int32), Errors.LineTooLong) # Read fields incrementally. `i` is the index to read from, beginning at 1. # The parse_x functions will return the next index to begin from (after the tab) (qname, i) = @? parse_str(mem, 1) (qlen, i) = @? parse_int(mem, i, false) (qstart, i) = @? parse_int(mem, i, true) - (qend, i) = @? parse_int(mem, i, false) + # qend may only be zero for unmapped records, e.g. minimap2's `--paf-no-hit` + # lines, so peek at the strand to check that before parsing it below. + (qend, i) = @? parse_int(mem, i, true) + if iszero(qend) && (i > lastindex(mem) || @inbounds(mem[i]) != UInt8('*')) + return ParserException((i - 1) % Int32, Errors.InvalidZero) + end qend ≥ qstart || return ParserException(i % Int32, Errors.BackwardsIndices) qend > qlen && return ParserException(i % Int32, Errors.PositionOutOfBounds) @@ -571,7 +578,9 @@ function parse_int( allow_zero::Bool, at_end::Bool = false, )::Union{Tuple{Int32, Int}, ParserException} - n = Int32(0) + # Accumulate in an Int, which can't overflow because we check against + # typemax(Int32) after every digit. + n = 0 i = from for outer i in from:lastindex(v) b = @inbounds v[i] @@ -580,11 +589,11 @@ function parse_int( i == from && return ParserException(i % Int32, Errors.EmptyInteger) (!allow_zero & iszero(n)) && return ParserException(i % Int32, Errors.InvalidZero) - return (n, i + 1) + return (n % Int32, i + 1) end b ∈ 0x30:0x39 || return ParserException(i % Int32, Errors.BadInteger) - n < 0 && return ParserException(i % Int32, Errors.IntegerOverflow) - n = Int32(10) * n + (b - Int32(48)) + n = 10 * n + (b - 0x30) + n > typemax(Int32) && return ParserException(i % Int32, Errors.IntegerOverflow) end # Note: It is not possible to get an InvalidZero error here, because # the last field (mapq) does allow zero, and for any other fields, @@ -596,7 +605,7 @@ function parse_int( elseif from > lastindex(v) ParserException(i % Int32, Errors.EmptyInteger) else - (n, i + 1) + (n % Int32, i + 1) end end @@ -780,18 +789,25 @@ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException parse_mem = buffer[1:last_index] # If there is no memory to parse (and we have already verified that # there are no more bytes in the IO), we end + # This is only the end of the file if the rest of the stream is whitespace. + # The buffer may end after the empty line, so keep filling it until EOF. if isempty(parse_mem) line = reader.line col = 0 - @inbounds for b in buffer - col += 1 - line += b == UInt8('\n') - col *= b != UInt8('\n') - if b ∉ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) - return ParserException(Errors.EmptyNonTailingLine, col % Int32, line) + n_scanned = 0 + while true + @inbounds for b in buffer[(n_scanned + 1):end] + col += 1 + line += b == UInt8('\n') + col *= b != UInt8('\n') + if b ∉ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) + return ParserException(Errors.EmptyNonTailingLine, col % Int32, line) + end end + n_scanned = length(buffer) + iszero(fill_buffer(io)) && return nothing + buffer = get_buffer(io) end - return nothing end # If this results in an error, we return the error. res = parse_line!(reader.record, parse_mem) diff --git a/test/runtests.jl b/test/runtests.jl index 38da693..17c7b4b 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -6,7 +6,24 @@ using Test using FormatSpecimens using MemoryViews using Aqua -using BufferIO: CursorReader +using BufferIO: CursorReader, AbstractBufReader +import BufferIO: get_buffer, fill_buffer, consume + +# A reader that delivers data in fixed chunks, like a pipe would +mutable struct ChunkReader <: AbstractBufReader + chunks::Vector{Vector{UInt8}} + buf::Vector{UInt8} +end +ChunkReader(chunks) = ChunkReader(map(c -> collect(codeunits(c)), chunks), UInt8[]) +get_buffer(r::ChunkReader) = ImmutableMemoryView(r.buf) +function fill_buffer(r::ChunkReader) + isempty(r.chunks) && return 0 + c = popfirst!(r.chunks) + append!(r.buf, c) + return length(c) +end +consume(r::ChunkReader, n::Int) = (deleteat!(r.buf, 1:n); nothing) +Base.close(::ChunkReader) = nothing using XAMAuxData: XAMAuxData using XAMAuxData.SAM: AuxTag @@ -114,6 +131,12 @@ end "rl" => 0, ) + # Unmapped record as written by minimap2 --paf-no-hit, with zeros + rec = parse(PAFRecord, "read1\t1523\t0\t0\t*\t*\t0\t0\t0\t0\t0\t0") + @test !is_mapped(rec) + @test rec.qname == "read1" + @test rec.qlen == 1523 + # Unmapped record with too few fieds @test try_parse("query\t123\t1\t100\t*\t*\t5\t19").kind == Errors.TooFewFields end @@ -170,6 +193,13 @@ end # Mapq too high @test with_replaced(good, 12, "256") == Errors.IntegerOverflow + + # Integer overflow, including on the last digit + @test with_replaced(good, 10, "2147483647").matches == typemax(Int32) + @test with_replaced(good, 10, "2147483648") == Errors.IntegerOverflow + @test with_replaced(good, 10, "9999999999") == Errors.IntegerOverflow + @test with_replaced(good, 12, "4294967296") == Errors.IntegerOverflow + @test with_replaced(good, 2, "2147483648") == Errors.IntegerOverflow end @testset "Bad reader" begin @@ -261,6 +291,16 @@ end @test first(reader, 2) isa Vector{PAFRecord} @test try_next!(reader).kind == Errors.EmptyNonTailingLine end + + # Blank line in the middle, where the buffer ends right after the blank line + good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" + reader = PAFReader(ChunkReader([good * "\n\n", good * "\n" * good * "\n"])) + @test try_next!(reader) isa PAFRecord + @test try_next!(reader).kind == Errors.EmptyNonTailingLine + + # Trailing whitespace delivered in multiple chunks is still fine + reader = PAFReader(ChunkReader([good * "\n\n", "\r\n", "\n"])) + @test length(collect(reader)) == 1 end @testset "Aqua" begin From 3c3127969b6215edb96e625714b66dbb589d15a6 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Sun, 18 Oct 2026 21:12:28 +0200 Subject: [PATCH 06/20] Make nonsense data access throw for unmapped records Trying to access e.g. `tlen` for an unmapped record now throws an `UnmappedError`. This simplifies some return types and streamlines how unmapped records handle garbage data. --- CHANGELOG.md | 5 ++ docs/src/index.md | 44 ++++++-------- src/PairwiseMappingFormat.jl | 115 +++++++++++++++++++++++++---------- test/runtests.jl | 12 +++- 4 files changed, 118 insertions(+), 58 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index b932232..700a069 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -27,6 +27,11 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 so the auxiliary data of records can no longer be modified. - **Breaking:** Reading from a closed `PAFReader` now throws, instead of behaving as if the reader is empty. +- **Breaking:** Accessing any property other than `qname` and `qlen` of an unmapped + record now throws the new exception `PairwiseMappingFormat.UnmappedError`, as does calling `query_coverage`, `target_coverage` + or `aln_identity` on it. Previously, `tname` and `is_rc` returned `nothing`, and the + remaining properties returned arbitrary values. As a result, `tname` is now always + a `StringView`, and `is_rc` always a `Bool`. ### Removed - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. diff --git a/docs/src/index.md b/docs/src/index.md index 8a7f6a0..968787f 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -89,8 +89,8 @@ underlying fields as they are stored in the `PAFRecord`. Names of the query sequences When the sequences are from a FASTA file, this is typically the _identifier_ of the FASTA records, i.e. the part of the header before the first space. -* `tname::Union{Nothing, StringView{A}}`. Same as `qname`, but for the target (i.e. subject) - sequence. May be `nothing` if the record is unmapped. +* `tname::StringView{A}`. Same as `qname`, but for the target (i.e. subject) + sequence. * `qlen` and `tlen`. Of type `Int`. Length of the full query and target sequence. * `qstart` and `tstart`. Of type `Int`. Leftmost (i.e. lowest) position of the alignment in the query and target, respectively. This does not take strandedness into account. @@ -104,7 +104,7 @@ underlying fields as they are stored in the `PAFRecord`. * `mapq::Union{Int, Nothing}`. Mapping quality, from 0:254, where 254 is the better quality. When calibrated, the probability the mapping is correct should be 10^(-mapq / 10). A value of `nothing` indicates the mapping quality is unavailable. -* `is_rc::Union{Bool, Nothing}` indicates whether the alignment matches the target on the forward (`false`) or reverse-complement `true` strand. A missing alignment is `nothing`. +* `is_rc::Bool` indicates whether the alignment matches the target on the forward (`false`) or reverse-complement `true` strand. The PAF format does not contain fields for alignment identity, query coverage or target coverage. @@ -124,41 +124,34 @@ record = parse(PAFRecord, data) ``` -Rarely, PAF records are unmapped. An unmapped record contain the query id, the query length, and any [auxiliary fields](@ref aux). -The target name and the strand of an unmapped return is `nothing`. -You can check if a record is mapped or unmapped with the function `is_mapped`: +Some PAF records are unmapped, for example those written by minimap2 with the +`--paf-no-hit` flag. An unmapped record only contains the query name, the query length, +and any [auxiliary fields](@ref aux). +You can check if a record is mapped or unmapped with the function [`is_mapped`](@ref): ```jldoctest unmapped -julia> unmapped = parse(PAFRecord, "myquery\t5032251\t9\t11\t*\t*\t5\t2\t1\t7\t4\t0"); +julia> unmapped = parse(PAFRecord, "myquery\t5032251\t0\t0\t*\t*\t0\t0\t0\t0\t0\t0"); -julia> @assert !is_mapped(unmapped) +julia> is_mapped(unmapped) +false julia> println(unmapped.qname) myquery julia> println(unmapped.qlen) 5032251 - -julia> @assert isnothing(unmapped.tname) - -julia> # The strand is given by is_rc, and is nothing for unmapped records - @assert isnothing(unmapped.is_rc) ``` -All other properties of unmapped records may contain arbitrary data, and should not be relied on. -For this reason, the functions [`target_coverage`](@ref), [`query_coverage`](@ref) and [`aln_identity`](@ref) -may give nonsense answers for unmapped records: +Accessing any other property of an unmapped record throws a [`PairwiseMappingFormat.UnmappedError`](@ref), +as does calling [`target_coverage`](@ref), [`query_coverage`](@ref) or [`aln_identity`](@ref) on it: ```jldoctest unmapped -julia> # Note! These results are arbitrary and not guaranteed to be stable - # across minor releases of this package - println(target_coverage(unmapped)) -Inf - -julia> println(query_coverage(unmapped)) -1.9871822768776836e-7 +julia> unmapped.tname +ERROR: UnmappedError: Cannot get `tname` of an unmapped PAFRecord. Check `is_mapped` first. +[...] -julia> println(aln_identity(unmapped)) -NaN +julia> query_coverage(unmapped) +ERROR: UnmappedError: Cannot get `query_coverage` of an unmapped PAFRecord. Check `is_mapped` first. +[...] ``` ## Low-level interface @@ -292,4 +285,5 @@ is_mapped PairwiseMappingFormat.Errors PairwiseMappingFormat.try_parse PairwiseMappingFormat.ParserException +PairwiseMappingFormat.UnmappedError ``` diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 32c3f5c..1c16ff3 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -9,7 +9,7 @@ export PAFReader, PAFRecord, aln_identity, aux_data, try_next!, query_coverage, target_coverage, is_mapped # Not exported, because this very PAF-specific, and operates on foreign types. -public try_parse, ParserException, Errors, Err +public try_parse, ParserException, UnmappedError, Errors, Err using PrecompileTools: @compile_workload, @setup_workload, @recompile_invalidations @@ -51,11 +51,12 @@ julia> record.qlen See also: [`PAFReader`](@ref) # Extended help -The following properties may be used: +The following properties may be used. For unmapped records, only `qname` and `qlen` +are available, and accessing any other property throws an [`UnmappedError`](@ref). +Use [`is_mapped`](@ref) to check if a record is mapped. * `qname::StringView`. The query name, May be empty, and contain any bytes except `\\t` and `\\n`. -* `tname::Union{StringView, Nothing}`. Target name. Like `qname`, but is `nothing` - if and only if the record is unmapped. +* `tname::StringView`. Target name. Like `qname`. * `qlen` and `tlen::Int`, and gives the length of the query and target sequences, respectively. This must be > 0. * `qstart`, `qend`, `tstart` and `tend::Int`, and give the starting and @@ -69,8 +70,8 @@ The following properties may be used: * `alnlen::Int` gives the length of the alignment * `mapq::Union{Int, Nothing}` gives the mapping quality, or `nothing` if this information is unavailable. -* `is_rc::Union{Bool,Nothing}` tells if the query and target strands are reverse-complement - relative to each other. Is `nothing` if the record is unmapped. +* `is_rc::Bool` tells if the query and target strands are reverse-complement + relative to each other. """ mutable struct PAFRecord # The line the record was parsed from, without the trailing newline. @@ -83,13 +84,15 @@ mutable struct PAFRecord qlen::Int32 # always >= 1 # Note: qstart/qend (and tstart/tend) are zero-indexed semi-open interval # in PAF format, but stored in this struct one-based. - qstart::Int32 # always >= 1, except when unmapped + # The fields below, and tname, are arbitrary when unmapped, and not accessible + # through the properties of the record. + qstart::Int32 # always >= 1 qend::Int32 # always >= qstart - tlen::Int32 # always >= 1, except when unmapped - tstart::Int32 # always >= 1, except when unmapped + tlen::Int32 # always >= 1 + tstart::Int32 # always >= 1 tend::Int32 # always >= tstart matches::Int32 - alnlen::Int32 # always >= 1, except when unmapped + alnlen::Int32 # always >= 1 mapq::UInt8 # stored as 0xff for missing strand::UInt8 # 0x00: *, 0x01: -, 0x02: + # 16 bits of padding, which we can use for something else later @@ -97,16 +100,24 @@ end # Since the interface of the records are defined in terms of its # properties, we need to be able to abstract over them -function Base.getproperty(record::PAFRecord, sym::Symbol) +# Note: This is inlined so `sym` is constant folded, leaving only the branch for +# that property, and the cheap check for mapped-ness. +@inline function Base.getproperty(record::PAFRecord, sym::Symbol) return if sym === :qname qname(record) - elseif sym == :tname + elseif sym === :qlen + Int(getfield(record, :qlen)) + elseif sym === :tname + check_mapped(record, sym) tname(record) - elseif sym == :mapq + elseif sym === :mapq + check_mapped(record, sym) mapq(record) - elseif sym in (:qlen, :qstart, :qend, :tlen, :tstart, :tend, :matches, :alnlen) + elseif sym in (:qstart, :qend, :tlen, :tstart, :tend, :matches, :alnlen) + check_mapped(record, sym) Int(getfield(record, sym)) - elseif sym == :is_rc + elseif sym === :is_rc + check_mapped(record, sym) is_rc(record) else error(lazy"Type PAFRecord has no property $sym") @@ -116,28 +127,68 @@ end """ is_mapped(x::PAFRecord) -> Bool -Compute whether the `PAFRecord` is mapped. An unmapped record will have -the properties `is_rc` and `tname` unavailable. -The properties `qname` and `qlen`, and the auxiliary data of an unmapped record -can be relied on, but the remaining properties contain arbitrary data. +Compute whether the `PAFRecord` is mapped. +Only the properties `qname` and `qlen`, and the auxiliary data of an unmapped record +are available. Accessing any other property of an unmapped record, or calling +[`query_coverage`](@ref), [`target_coverage`](@ref) or [`aln_identity`](@ref) on it, +throws an [`UnmappedError`](@ref). Note that some PAF parsers do not handle unmapped records correctly, so be wary when writing unmapped records. # Examples -```jldoctes +```jldoctest julia> is_mapped(record) true + +julia> unmapped = parse(PAFRecord, "read1\\t1523\\t0\\t0\\t*\\t*\\t0\\t0\\t0\\t0\\t0\\t0"); + +julia> is_mapped(unmapped) +false ``` See also: [`PAFRecord`](@ref) """ is_mapped(x::PAFRecord) = !iszero(getfield(x, :strand)) -function is_rc(record::PAFRecord)::Union{Bool, Nothing} - x = getfield(record, :strand) - return iszero(x) ? nothing : x == 0x02 +""" + UnmappedError + +The exception thrown when accessing a property of an unmapped [`PAFRecord`](@ref) which +only exists for mapped records, such as `tname`, or when calling a function which +requires a mapped record, such as [`query_coverage`](@ref). +The name of the property or function may be accessed with the `.name` field. + +# Examples +```jldoctest +julia> unmapped = parse(PAFRecord, "read1\\t1523\\t0\\t0\\t*\\t*\\t0\\t0\\t0\\t0\\t0\\t0"); + +julia> err = try unmapped.tname catch e e end; + +julia> err isa PAF.UnmappedError +true + +julia> err.name +:tname +``` + +See also: [`is_mapped`](@ref) +""" +struct UnmappedError <: Exception + name::Symbol +end + +function Base.showerror(io::IO, err::UnmappedError) + return print(io, "UnmappedError: Cannot get `", err.name, "` of an unmapped PAFRecord. Check `is_mapped` first.") end +# Keep the error path out of line, so the check is cheap to inline +@noinline throw_unmapped(name::Symbol) = throw(UnmappedError(name)) + +check_mapped(record::PAFRecord, name::Symbol) = (is_mapped(record) || throw_unmapped(name); nothing) + +# Note: The functions below do not check if the record is mapped. +is_rc(record::PAFRecord) = getfield(record, :strand) == 0x02 + # For tab completion function Base.propertynames(::PAFRecord) return ( @@ -198,14 +249,10 @@ function qname(record::PAFRecord) end # Return StringView to not allocate -function tname(record::PAFRecord)::Union{StringView, Nothing} - return if is_mapped(record) - start = getfield(record, :tname_start) - span = start:(start + getfield(record, :tname_len) - 1) - StringView(getfield(record, :data)[span]) - else - nothing - end +function tname(record::PAFRecord) + start = getfield(record, :tname_start) + span = start:(start + getfield(record, :tname_len) - 1) + return StringView(getfield(record, :data)[span]) end # Note: We return nothing instead of missing because fuck missing @@ -220,6 +267,7 @@ end Return the nucleotide / residue identity, defined as the number of matches divided by the alignment length. This is a number between 0 and 1. +Throws an [`UnmappedError`](@ref) if the record is unmapped. # Examples ```jldoctest @@ -230,6 +278,7 @@ julia> aln_identity(record) See also: [`query_coverage`](@ref), [`target_coverage`](@ref) """ function aln_identity(record::PAFRecord)::Float64 + check_mapped(record, :aln_identity) return record.matches / record.alnlen end @@ -237,6 +286,7 @@ end query_coverage(rec::Record) -> Float64 Compute the fraction of the query covered by the alignment. +Throws an [`UnmappedError`](@ref) if the record is unmapped. # Examples ```jldoctest @@ -247,6 +297,7 @@ julia> query_coverage(record) See also: [`target_coverage`](@ref), [`aln_identity`](@ref) """ function query_coverage(record::PAFRecord)::Float64 + check_mapped(record, :query_coverage) return (record.qend - record.qstart + Int32(1)) / record.qlen end @@ -254,6 +305,7 @@ end target_coverage(rec::Record) -> Float64 Compute the fraction of the target covered by the alignment. +Throws an [`UnmappedError`](@ref) if the record is unmapped. # Examples ```jldoctest @@ -264,6 +316,7 @@ julia> target_coverage(record) See also: [`query_coverage`](@ref), [`aln_identity`](@ref) """ function target_coverage(record::PAFRecord)::Float64 + check_mapped(record, :target_coverage) return (record.tend - record.tstart + Int32(1)) / record.tlen end diff --git a/test/runtests.jl b/test/runtests.jl index 17c7b4b..13c76ed 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -63,6 +63,10 @@ end @test rec1.alnlen == 301142 @test rec1.mapq == 0 + # tname and is_rc are type stable, since they throw for unmapped records + @test (@inferred (r -> r.tname)(rec1)) isa AbstractString + @test (@inferred (r -> r.is_rc)(rec1)) isa Bool + @test_throws Exception rec1.foobar @test_throws Exception rec1.strand @test_throws Exception rec1.data @@ -118,9 +122,13 @@ end # Empty record rec = PAFReader(first, open(joinpath(path_of_format("PAF"), "good3.paf"))) - @test isnothing(rec.tname) - @test isnothing(rec.is_rc) @test !is_mapped(rec) + for sym in (:tname, :is_rc, :mapq, :qstart, :qend, :tlen, :tstart, :tend, :matches, :alnlen) + @test_throws PairwiseMappingFormat.UnmappedError(sym) getproperty(rec, sym) + end + for f in (query_coverage, target_coverage, aln_identity) + @test_throws PairwiseMappingFormat.UnmappedError(Symbol(f)) f(rec) + end @test rec.qname == "query" @test rec.qlen == 56 @test Dict(aux_data(rec)) == Dict{AuxTag, Any}( From 9d711bb79e1e9f91a5a426113ba54d2e34ec1017 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Sun, 18 Oct 2026 21:18:15 +0200 Subject: [PATCH 07/20] Switch internal fields to UInt32 Doubles the length of accepted templates for very little cost --- CHANGELOG.md | 2 ++ docs/src/index.md | 1 + src/PairwiseMappingFormat.jl | 69 +++++++++++++++++++----------------- test/runtests.jl | 19 ++++++++-- 4 files changed, 56 insertions(+), 35 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 700a069..a5fe346 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -32,6 +32,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 or `aln_identity` on it. Previously, `tname` and `is_rc` returned `nothing`, and the remaining properties returned arbitrary values. As a result, `tname` is now always a `StringView`, and `is_rc` always a `Bool`. +- Sequence lengths and positions up to `typemax(UInt32)` (about 4.3 Gbp) are now supported, + up from `typemax(Int32)`. Larger integers are rejected with `IntegerOverflow`. ### Removed - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. diff --git a/docs/src/index.md b/docs/src/index.md index 968787f..358e0a8 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -92,6 +92,7 @@ underlying fields as they are stored in the `PAFRecord`. * `tname::StringView{A}`. Same as `qname`, but for the target (i.e. subject) sequence. * `qlen` and `tlen`. Of type `Int`. Length of the full query and target sequence. + Lengths and positions up to `typemax(UInt32)`, about 4.3 Gbp, are supported. * `qstart` and `tstart`. Of type `Int`. Leftmost (i.e. lowest) position of the alignment in the query and target, respectively. This does not take strandedness into account. * `qend` and `tend`. Of type `Int`. Rightmost (i.e. highest) position of the alignment diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 1c16ff3..6e4d3da 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -58,7 +58,7 @@ Use [`is_mapped`](@ref) to check if a record is mapped. except `\\t` and `\\n`. * `tname::StringView`. Target name. Like `qname`. * `qlen` and `tlen::Int`, and gives the length of the query and target - sequences, respectively. This must be > 0. + sequences, respectively. This must be > 0, and at most `typemax(UInt32)` (about 4.3 Gbp). * `qstart`, `qend`, `tstart` and `tend::Int`, and give the starting and ending positions of the alignments on the query and target strand. These uses one-based, closed interval indexing as usual in Julia, but unlike @@ -81,18 +81,21 @@ mutable struct PAFRecord tname_start::Int32 tname_len::Int32 # always zero when unmapped aux_start::Int32 # lastindex(data) + 1 if there is no aux data - qlen::Int32 # always >= 1 - # Note: qstart/qend (and tstart/tend) are zero-indexed semi-open interval - # in PAF format, but stored in this struct one-based. + # Positions and lengths are stored as UInt32, so sequences of up to + # typemax(UInt32) bytes are supported. + qlen::UInt32 # always >= 1 + # Note: qstart/qend (and tstart/tend) are zero-indexed semi-open intervals + # in the PAF format, and are stored as such. The properties are one-based. + # The start is stored zero-based, such that the conversion can't overflow. # The fields below, and tname, are arbitrary when unmapped, and not accessible # through the properties of the record. - qstart::Int32 # always >= 1 - qend::Int32 # always >= qstart - tlen::Int32 # always >= 1 - tstart::Int32 # always >= 1 - tend::Int32 # always >= tstart - matches::Int32 - alnlen::Int32 # always >= 1 + qstart::UInt32 + qend::UInt32 # always >= qstart + tlen::UInt32 # always >= 1 + tstart::UInt32 + tend::UInt32 # always >= tstart + matches::UInt32 + alnlen::UInt32 # always >= 1 mapq::UInt8 # stored as 0xff for missing strand::UInt8 # 0x00: *, 0x01: -, 0x02: + # 16 bits of padding, which we can use for something else later @@ -113,7 +116,10 @@ end elseif sym === :mapq check_mapped(record, sym) mapq(record) - elseif sym in (:qstart, :qend, :tlen, :tstart, :tend, :matches, :alnlen) + elseif sym in (:qstart, :tstart) + check_mapped(record, sym) + Int(getfield(record, sym)) + 1 + elseif sym in (:qend, :tlen, :tend, :matches, :alnlen) check_mapped(record, sym) Int(getfield(record, sym)) elseif sym === :is_rc @@ -209,7 +215,7 @@ end function PAFRecord() data = ImmutableMemoryView(Memory{UInt8}()) - return PAFRecord(data, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0xff, 0x00) + return PAFRecord(data, 0, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 0, 0xff, 0x00) end function Base.show(io::IO, ::MIME"text/plain", record::PAFRecord) @@ -298,7 +304,7 @@ See also: [`target_coverage`](@ref), [`aln_identity`](@ref) """ function query_coverage(record::PAFRecord)::Float64 check_mapped(record, :query_coverage) - return (record.qend - record.qstart + Int32(1)) / record.qlen + return (record.qend - record.qstart + 1) / record.qlen end """ @@ -317,7 +323,7 @@ See also: [`query_coverage`](@ref), [`aln_identity`](@ref) """ function target_coverage(record::PAFRecord)::Float64 check_mapped(record, :target_coverage) - return (record.tend - record.tstart + Int32(1)) / record.tlen + return (record.tend - record.tstart + 1) / record.tlen end """ @@ -548,21 +554,20 @@ function parse_line!( end # Fill in fields of the struct. The aux data begins at `i`. - # Note: PAF uses zero-based semi-open indexing like e.g. Python, - # so compensate for that in qstart and qend. + # The integers have been checked to fit in an UInt32 by `parse_int`. record.data = mem record.qname_len = length(qname) % Int32 record.tname_start = first(tname) % Int32 record.tname_len = length(tname) % Int32 record.aux_start = i % Int32 - record.qlen = qlen - record.qstart = qstart + one(qstart) - record.qend = qend - record.tlen = tlen - record.tstart = tstart + one(tstart) - record.tend = tend - record.matches = matches - record.alnlen = alnlen + record.qlen = qlen % UInt32 + record.qstart = qstart % UInt32 + record.qend = qend % UInt32 + record.tlen = tlen % UInt32 + record.tstart = tstart % UInt32 + record.tend = tend % UInt32 + record.matches = matches % UInt32 + record.alnlen = alnlen % UInt32 record.mapq = mapq record.strand = strand @@ -576,7 +581,7 @@ function finish_unmapped!( record::PAFRecord, mem::ImmutableMemoryView{UInt8}, qname::UnitRange{Int}, - qlen::Int32, + qlen::Int, i::Int, ) # Determine the position of aux data @@ -597,7 +602,7 @@ function finish_unmapped!( record.tname_start = 1 record.tname_len = 0 record.aux_start = aux_start % Int32 - record.qlen = qlen + record.qlen = qlen % UInt32 record.qstart = 0 record.qend = 0 record.tlen = 0 @@ -630,9 +635,9 @@ function parse_int( from::Int, allow_zero::Bool, at_end::Bool = false, - )::Union{Tuple{Int32, Int}, ParserException} + )::Union{Tuple{Int, Int}, ParserException} # Accumulate in an Int, which can't overflow because we check against - # typemax(Int32) after every digit. + # typemax(UInt32) after every digit. The result therefore fits in an UInt32. n = 0 i = from for outer i in from:lastindex(v) @@ -642,11 +647,11 @@ function parse_int( i == from && return ParserException(i % Int32, Errors.EmptyInteger) (!allow_zero & iszero(n)) && return ParserException(i % Int32, Errors.InvalidZero) - return (n % Int32, i + 1) + return (n, i + 1) end b ∈ 0x30:0x39 || return ParserException(i % Int32, Errors.BadInteger) n = 10 * n + (b - 0x30) - n > typemax(Int32) && return ParserException(i % Int32, Errors.IntegerOverflow) + n > typemax(UInt32) && return ParserException(i % Int32, Errors.IntegerOverflow) end # Note: It is not possible to get an InvalidZero error here, because # the last field (mapq) does allow zero, and for any other fields, @@ -658,7 +663,7 @@ function parse_int( elseif from > lastindex(v) ParserException(i % Int32, Errors.EmptyInteger) else - (n % Int32, i + 1) + (n, i + 1) end end diff --git a/test/runtests.jl b/test/runtests.jl index 13c76ed..746c387 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -203,11 +203,24 @@ end @test with_replaced(good, 12, "256") == Errors.IntegerOverflow # Integer overflow, including on the last digit - @test with_replaced(good, 10, "2147483647").matches == typemax(Int32) - @test with_replaced(good, 10, "2147483648") == Errors.IntegerOverflow + @test with_replaced(good, 10, "4294967295").matches == typemax(UInt32) + @test with_replaced(good, 10, "4294967296") == Errors.IntegerOverflow @test with_replaced(good, 10, "9999999999") == Errors.IntegerOverflow + @test with_replaced(good, 10, "99999999999999999999") == Errors.IntegerOverflow @test with_replaced(good, 12, "4294967296") == Errors.IntegerOverflow - @test with_replaced(good, 2, "2147483648") == Errors.IntegerOverflow + @test with_replaced(good, 2, "4294967296") == Errors.IntegerOverflow + + # Sequences longer than typemax(Int32), up to typemax(UInt32) + m = typemax(UInt32) + rec = try_parse("q\t$m\t$(m - 1)\t$m\t+\tt\t$m\t0\t$m\t1\t1\t0") + @test (rec.qlen, rec.qstart, rec.qend) == (m, m, m) + @test (rec.tlen, rec.tstart, rec.tend) == (m, 1, m) + @test query_coverage(rec) == 1 / m + @test target_coverage(rec) == 1.0 + + # Empty interval at the very end of the sequence doesn't overflow + rec = try_parse("q\t$m\t$m\t$m\t+\tt\t5\t0\t5\t1\t1\t0") + @test (rec.qstart, rec.qend) == (Int(m) + 1, m) end @testset "Bad reader" begin From 425b4a8de60f4e965a720a5174e70360349a737d Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Mon, 19 Oct 2026 14:55:38 +0200 Subject: [PATCH 08/20] Make `is_empty` not consume a record Also refactor checking empty file tail. --- AGENTS.md | 14 +++++ CHANGELOG.md | 5 ++ src/PairwiseMappingFormat.jl | 118 +++++++++++++++++++++-------------- test/runtests.jl | 46 ++++++++++++++ 4 files changed, 137 insertions(+), 46 deletions(-) create mode 100644 AGENTS.md diff --git a/AGENTS.md b/AGENTS.md new file mode 100644 index 0000000..fc2bbc7 --- /dev/null +++ b/AGENTS.md @@ -0,0 +1,14 @@ +# AGENTS.md +This file contains instructions for AI agents working with this repository. + +Purpose: PairwiseMappingFormat.jl provides parsing for the bioinformatics PAF format. + +Design goals: +* Be correct, even in edge cases. Unrealistic edge cases like 2 GiB lines or 4 GiB references may not be optimized for, but should be handled. +* Be performant + +Scope: +* Writing or mutating records is intentionally not supported yet. + +Guidelines for developing +* After making changes, format with `runic -i .`. If Runic is not installed, alert the user and do not attempt to install it. diff --git a/CHANGELOG.md b/CHANGELOG.md index a5fe346..53b965c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -37,3 +37,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Removed - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. + +### Fixed +- `isempty(reader)` no longer consumes a record from a `PAFReader`. `PAFReader` now + implements `Base.isdone`, which checks for more records without advancing the reader + or invalidating previously returned records. diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 6e4d3da..bc79df9 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -814,12 +814,34 @@ true ``` """ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} - io = reader.io - # Consume the previous line, now that the previous record is invalidated. + consume_previous!(reader) + (buffer, last_index, line_len) = buffer_line(reader.io) + # If the line is empty, this is the end of the stream, if the rest of it is whitespace. + iszero(last_index) && return check_tail(reader, buffer) + # If this results in an error, we return the error. + res = parse_line!(reader.record, buffer[1:last_index]) + res isa ParserException && + return ParserException(res.kind, getfield(res, :index_in_line), reader.line) + # Else, we update the reader. This means `next!` will continuously return + # an error if an error is found. + reader.to_consume = line_len + reader.line += 1 + return res +end + +# Consume the previous line. This invalidates the record if it is a view into the buffer. +function consume_previous!(reader::PAFReader) if !iszero(reader.to_consume) - consume(io, reader.to_consume) + consume(reader.io, reader.to_consume) reader.to_consume = 0 end + return nothing +end + +# Fill the buffer until it contains a full line, or EOF is reached. +# Return the buffer, the index of the last byte of the line before the \n (or \r\n), +# and the number of bytes in the line including the newline. +function buffer_line(io::ExpandReader)::Tuple{ImmutableMemoryView{UInt8}, Int, Int} buffer = get_buffer(io) newline = findfirst(==(0x0a), buffer) # Keep reading until we find a newline or reach EOF @@ -828,56 +850,60 @@ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException n_scanned = length(buffer) # This signals EOF - we might still have data in the buffer, if the underlying stream # ends with a record that doesn't end in a newline - iszero(fill_buffer(io)) && break + iszero(fill_buffer(io)) && return (buffer, length(buffer), length(buffer)) buffer = get_buffer(io) newline = findnext(==(0x0a), buffer, n_scanned + 1) end - has_windows_newline = false # Last index before \n, and before \r if the newline is \r\n as on Windows - last_index = if newline === nothing - length(buffer) - else - if isone(newline) - 0 - else - has_windows_newline = @inbounds(buffer[newline - 1]) == UInt8('\r') - newline - (1 + has_windows_newline) - end - end - parse_mem = buffer[1:last_index] - # If there is no memory to parse (and we have already verified that - # there are no more bytes in the IO), we end - # This is only the end of the file if the rest of the stream is whitespace. - # The buffer may end after the empty line, so keep filling it until EOF. - if isempty(parse_mem) - line = reader.line - col = 0 - n_scanned = 0 - while true - @inbounds for b in buffer[(n_scanned + 1):end] - col += 1 - line += b == UInt8('\n') - col *= b != UInt8('\n') - if b ∉ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) - return ParserException(Errors.EmptyNonTailingLine, col % Int32, line) - end + has_windows_newline = !isone(newline) && @inbounds(buffer[newline - 1]) == UInt8('\r') + return (buffer, newline - (1 + has_windows_newline), newline) +end + +# Called when the next line is empty. Return `nothing` if the rest of the stream is +# whitespace, else an error. +# The buffer may end after the empty line, so keep filling it until EOF. +function check_tail(reader::PAFReader, buffer::ImmutableMemoryView{UInt8})::Union{Nothing, ParserException} + io = reader.io + line = reader.line + col = 0 + n_scanned = 0 + while true + @inbounds for b in buffer[(n_scanned + 1):end] + col += 1 + line += b == UInt8('\n') + col *= b != UInt8('\n') + if b ∉ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) + return ParserException(Errors.EmptyNonTailingLine, col % Int32, line) end - n_scanned = length(buffer) - iszero(fill_buffer(io)) && return nothing - buffer = get_buffer(io) end + n_scanned = length(buffer) + iszero(fill_buffer(io)) && return nothing + buffer = get_buffer(io) end - # If this results in an error, we return the error. - res = parse_line!(reader.record, parse_mem) - res isa ParserException && - return ParserException(res.kind, getfield(res, :index_in_line), reader.line) - # Else, we update the reader. This means `next!` will continuously return - # an error if an error is found. - # If the newline was found, we compensate for the newline by also consuming that byte, - # hence the extra `newline !== nothing`. - reader.to_consume = last_index + has_windows_newline + (newline !== nothing) - reader.line += 1 - return res + return +end + +# This is also used by `isempty`. It does not advance the reader, and does not +# invalidate any records returned previously. +function Base.isdone(reader::PAFReader, _...) + # Fast path: If the next line begins with any byte other than a newline, it is not + # empty, so `try_next!` would return a record or an error. + buffer = get_buffer(reader.io) + n = reader.to_consume + if length(buffer) > n + b = @inbounds buffer[n + 1] + b ∈ (UInt8('\r'), UInt8('\n')) || return false + end + # Slow path: We need to consume the previous line, and fill the buffer, which may + # overwrite the bytes the reader's record is a view into. So, make the record own its + # data first. This keeps any references to the record held by users valid. + if !iszero(n) + record = reader.record + setfield!(record, :data, copy(getfield(record, :data))) + consume_previous!(reader) + end + (buffer, last_index, _) = buffer_line(reader.io) + return iszero(last_index) && check_tail(reader, buffer) === nothing end @setup_workload begin diff --git a/test/runtests.jl b/test/runtests.jl index 746c387..d86635f 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -324,6 +324,52 @@ end @test length(collect(reader)) == 1 end +@testset "isempty" begin + good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" + + # isempty does not consume records + reader = PAFReader(IOBuffer(join(fill(good, 3), '\n'))) + @test !isempty(reader) + @test !isempty(reader) + @test length(collect(reader)) == 3 + @test isempty(reader) + + # isempty does not invalidate the record with copy=false, even when it needs to + # fill the buffer + for chunks in ([good * "\n", good * "\n"], [good * "\n", "\n", "\r\n"], [good]) + reader = PAFReader(ChunkReader(chunks); copy = false) + rec = try_next!(reader) + expected = copy(rec) + while !isempty(reader) + @test rec.qname == expected.qname + @test rec.tname == expected.tname + @test aux_data(rec) == aux_data(expected) + rec = try_next!(reader) + @test rec isa PAFRecord + end + @test rec.qname == expected.qname + @test try_next!(reader) === nothing + end + + # isempty agrees with try_next! for all inputs and chunkings + inputs = [ + "", "\n", "\r\n", "\n\n \t\r\n", good, good * "\n", good * "\r\n\n", + good * "\n\n" * good, good * "\nbad", good * "\n \n", good * "\n\r", + " ", "\nV", good * "\n\n\n\nx", + ] + for data in inputs, split_at in 0:length(data) + # Empty chunks are not allowed, since fill_buffer returning zero signals EOF + chunks = filter!(!isempty, [data[1:split_at], data[(split_at + 1):end]]) + reader = PAFReader(ChunkReader(chunks)) + while true + done = isempty(reader) + res = try_next!(reader) + @test done == (res === nothing) + res isa PAFRecord || break + end + end +end + @testset "Aqua" begin Aqua.test_all(PairwiseMappingFormat; ambiguities = false) end From d94afdab11336b592dcb4982af9487a0dc6ce271 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Mon, 19 Oct 2026 15:04:20 +0200 Subject: [PATCH 09/20] Validate unmapped records --- CHANGELOG.md | 2 ++ src/PairwiseMappingFormat.jl | 35 ++++++++++++++++--------------- test/runtests.jl | 40 ++++++++++++++++++++++++++++++------ 3 files changed, 54 insertions(+), 23 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 53b965c..f353094 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -34,6 +34,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 a `StringView`, and `is_rc` always a `Bool`. - Sequence lengths and positions up to `typemax(UInt32)` (about 4.3 Gbp) are now supported, up from `typemax(Int32)`. Larger integers are rejected with `IntegerOverflow`. +- **Breaking:** Fields 6-12 of unmapped records are now validated like those of mapped + records, except that the integer fields may be zero. Previously, their content was ignored. ### Removed - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index bc79df9..cb72eb3 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -524,7 +524,7 @@ function parse_line!( return ParserException(lastindex(mem) % Int32, Errors.TooFewFields) b = @inbounds mem[i] strand = if b == UInt8('*') - return finish_unmapped!(record, mem, qname, qlen, i + 2) + 0x00 elseif b == UInt8('-') 0x02 elseif b == UInt8('+') @@ -534,6 +534,7 @@ function parse_line!( end @inbounds mem[i + 1] == UInt8('\t') || return ParserException((i + 1) % Int32, Errors.TooFewFields) + iszero(strand) && return finish_unmapped!(record, mem, qname, qlen, i + 2) # Load rest of the fields (tname, i) = @? parse_str(mem, i + 2) # i + 2 because strand plus tab took two bytes @@ -574,8 +575,9 @@ function parse_line!( return record end -# If the record is unmapped we don't even bother reading -# the rest of the fields, so take this fast path +# Parse the remaining fields of an unmapped record. +# These are validated like for mapped records, except that the integers may be zero, +# as written by e.g. minimap2's `--paf-no-hit`. Their values are not stored. # i is the byte index of the strand plus two function finish_unmapped!( record::PAFRecord, @@ -583,25 +585,24 @@ function finish_unmapped!( qname::UnitRange{Int}, qlen::Int, i::Int, - ) - # Determine the position of aux data - for _ in 1:6 - i = findnext(==(UInt8('\t')), mem, i) - i = if isnothing(i) - return ParserException(lastindex(mem) % Int32, Errors.TooFewFields) - else - i + 1 - end - end - i = findnext(==(UInt8('\t')), mem, i) - aux_start = isnothing(i) ? lastindex(mem) + 1 : i + 1 + )::Union{PAFRecord, ParserException} + (_, i) = @? parse_str(mem, i) + (tlen, i) = @? parse_int(mem, i, true) + (tstart, i) = @? parse_int(mem, i, true) + (tend, i) = @? parse_int(mem, i, true) + tend ≥ tstart || return ParserException(i % Int32, Errors.BackwardsIndices) + tend > tlen && return ParserException(i % Int32, Errors.PositionOutOfBounds) + (_, i) = @? parse_int(mem, i, true) # matches + (_, i) = @? parse_int(mem, i, true) # alnlen + (mapq, i) = @? parse_int(mem, i, true, true) + mapq > 255 && return ParserException(i % Int32, Errors.IntegerOverflow) - # Fill in fields + # Fill in fields. The aux data begins at `i`. record.data = mem record.qname_len = length(qname) % Int32 record.tname_start = 1 record.tname_len = 0 - record.aux_start = aux_start % Int32 + record.aux_start = i % Int32 record.qlen = qlen % UInt32 record.qstart = 0 record.qend = 0 diff --git a/test/runtests.jl b/test/runtests.jl index d86635f..d6fab10 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -106,6 +106,17 @@ end @test isapprox(target_coverage(get_rec(17, 16, 15, 16)), 6 / 16) end +function with_replaced_str(s::String, field::Integer, new::String) + fields = split(s, '\t') + fields[field] = new + return join(fields, "\t") +end + +function with_replaced(s::String, field::Integer, new::String) + y = try_parse(with_replaced_str(s, field, new)) + return y isa PairwiseMappingFormat.ParserException ? y.kind : y +end + @testset "Edge case records" begin (rec1, rec2, rec3) = PAFReader(open(joinpath(path_of_format("PAF"), "good2.paf"))) do reader y = first(reader, 3) @@ -147,13 +158,30 @@ end # Unmapped record with too few fieds @test try_parse("query\t123\t1\t100\t*\t*\t5\t19").kind == Errors.TooFewFields -end -function with_replaced(s::String, field::Integer, new::String) - fields = split(s, '\t') - fields[field] = new - y = try_parse(join(fields, "\t")) - return y isa PairwiseMappingFormat.ParserException ? y.kind : y + # The fields of unmapped records are validated like for mapped records, + # except that integers may be zero + unmapped = "read1\t1523\t0\t0\t*\t*\t0\t0\t0\t0\t0\t0" + for field in 6:12 + @test with_replaced(unmapped, field, "0") isa PAFRecord + end + # Use a longer target name, so emptying a field does not make the line too short + long_unmapped = with_replaced_str(unmapped, 6, "some_target_name") + for field in 7:12 + @test with_replaced(long_unmapped, field, "") == Errors.EmptyInteger + @test with_replaced(unmapped, field, "-1") == Errors.BadInteger + @test with_replaced(unmapped, field, "abc") == Errors.BadInteger + @test with_replaced(unmapped, field, "4294967296") == Errors.IntegerOverflow + end + @test with_replaced(unmapped, 12, "255") isa PAFRecord + @test with_replaced(unmapped, 12, "256") == Errors.IntegerOverflow + @test with_replaced(unmapped, 9, "1") == Errors.PositionOutOfBounds + @test try_parse("read1\t1523\t0\t0\t*\t*\t10\t5\t4\t0\t0\t0").kind == Errors.BackwardsIndices + @test try_parse("read1\t1523\t0\t0\t*\t*\t10\t0\t5\t0\t0\t0") isa PAFRecord + @test try_parse("read1\t1523\t0\t0\t*x*\t10\t0\t5\t0\t0\t0").kind == Errors.TooFewFields + # Aux data is found after the mapq field + rec = try_parse("read1\t1523\t0\t0\t*\t*\t0\t0\t0\t0\t0\t0\tcm:i:4") + @test Dict(aux_data(rec)) == Dict{AuxTag, Any}("cm" => 4) end @testset "Errors" begin From bfec1e811a4787d98c001c25757284ac50cbfe89 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Mon, 19 Oct 2026 15:37:51 +0200 Subject: [PATCH 10/20] Avoid invalidating records on `isdone` --- CHANGELOG.md | 12 ++- docs/src/index.md | 26 ++++-- src/PairwiseMappingFormat.jl | 176 ++++++++++++++++++++--------------- test/runtests.jl | 68 +++++++++++--- 4 files changed, 178 insertions(+), 104 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index f353094..4cb322f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -23,6 +23,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - `try_parse` and `parse` return a view into the input, so mutating the input also mutates the record. - `copy(record)` returns a record that owns its data. +- **Breaking:** `PAFRecord` is now immutable. With `copy=false`, iterating a `PAFReader` + no longer returns the same object on each iteration, but a new record each time. + On Julia versions before 1.14, returning a `PAFRecord` from a non-inlined function, + such as `try_next!` or `iterate`, allocates. - **Breaking:** `aux_data` now returns an immutable `Auxiliary{ImmutableMemoryView{UInt8}}`, so the auxiliary data of records can no longer be modified. - **Breaking:** Reading from a closed `PAFReader` now throws, instead of @@ -41,6 +45,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. ### Fixed -- `isempty(reader)` no longer consumes a record from a `PAFReader`. `PAFReader` now - implements `Base.isdone`, which checks for more records without advancing the reader - or invalidating previously returned records. +- `isempty(reader)` no longer consumes a record from a `PAFReader`. It may invalidate + the most recently returned record if that record is a view into the reader's buffer, + i.e. if it was returned by `try_next!`, or iterated with `copy=false`. + `PAFReader` now implements `Base.isdone`, which never invalidates records, and returns + `missing` when it can't determine whether the reader is done without doing so. diff --git a/docs/src/index.md b/docs/src/index.md index 358e0a8..70fcb4c 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -47,24 +47,29 @@ julia> PAFReader(open(path_to_paf)) do reader ``` The [`PAFReader`](@ref) constructor takes an optional keyword `copy`, which defaults to `true`. -If it is `false`, the record will iterate the _same_ [`PAFRecord`](@ref) object, overwriting -it in-place. This record is a view into the reader's buffer, which is only valid until the -reader is advanced. +If it is `false`, the iterated [`PAFRecord`](@ref)s are views into the reader's buffer, +instead of owning a copy of their data. +Such a record is only valid until the reader is advanced, and may be invalidated by +calling [`isempty`](@ref Base.isempty(::PAFReader)) on the reader. This reduces allocations and copying and gives a slight speed increase, but may result -in bugs if records, or references of records of previous iterations are stored and -unexpectedly overwritten: +in bugs if records from previous iterations are stored, since their data may be +silently overwritten. Use `copy` on a record to obtain a record that owns its data: ```jldoctest -julia> records = PAFReader(collect, open(path_to_paf); copy=false); - -julia> println(map(i -> i.qlen, records)) # NB: All the same record! -[288659, 288659, 288659] +julia> PAFReader(open(path_to_paf); copy=false) do reader + [copy(record) for record in reader if record.qlen > 290000] + end |> length +2 ``` +!!! note + On Julia versions before 1.14, returning a `PAFRecord` from a function that is not + inlined allocates, even with `copy=false`. + At the moments, readers do not support seeking. ## Record -The mutable [`PAFRecord`](@ref) object represents a single line in a PAF file. +The immutable [`PAFRecord`](@ref) object represents a single line in a PAF file. The individual columns of the PAF line is obtained by accessing properties of the records: ```jldoctest record @@ -282,6 +287,7 @@ aln_identity query_coverage target_coverage try_next! +Base.isempty(::PAFReader) is_mapped PairwiseMappingFormat.Errors PairwiseMappingFormat.try_parse diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index cb72eb3..5206548 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -30,7 +30,7 @@ using BufferIO: AbstractBufReader, BufReader, ExpandReader, get_buffer, fill_buf """ PAFRecord -Mutable type representing a single PAF line. The content of the record +Immutable type representing a single PAF line. The content of the record is accessed through its properties. A `PAFRecord` does not necessarily own its data, but may be a view into the @@ -73,7 +73,7 @@ Use [`is_mapped`](@ref) to check if a record is mapped. * `is_rc::Bool` tells if the query and target strands are reverse-complement relative to each other. """ -mutable struct PAFRecord +struct PAFRecord # The line the record was parsed from, without the trailing newline. # May be a view into memory not owned by the record, e.g. a reader's buffer. data::ImmutableMemoryView{UInt8} @@ -490,17 +490,14 @@ PairwiseMappingFormat.ParserException See also: [`PAFRecord`](@ref), [`try_next!`](@ref) """ try_parse(x) = try_parse(ImmutableMemoryView(x)) -try_parse(x::ImmutableMemoryView{UInt8}) = parse_line!(PAFRecord(), x) +try_parse(x::ImmutableMemoryView{UInt8}) = parse_line(x) function Base.parse(::Type{PAFRecord}, s::Union{String, SubString{String}}) y = try_parse(s) return y isa PAFRecord ? y : throw(y) end -function parse_line!( - record::PAFRecord, - mem::ImmutableMemoryView{UInt8}, - )::Union{PAFRecord, ParserException} +function parse_line(mem::ImmutableMemoryView{UInt8})::Union{PAFRecord, ParserException} # Note: The aux data may begin at lastindex(mem) + 1, which must fit in an Int32 lastindex(mem) >= typemax(Int32) && return ParserException(typemax(Int32), Errors.LineTooLong) @@ -534,7 +531,7 @@ function parse_line!( end @inbounds mem[i + 1] == UInt8('\t') || return ParserException((i + 1) % Int32, Errors.TooFewFields) - iszero(strand) && return finish_unmapped!(record, mem, qname, qlen, i + 2) + iszero(strand) && return finish_unmapped(mem, qname, qlen, i + 2) # Load rest of the fields (tname, i) = @? parse_str(mem, i + 2) # i + 2 because strand plus tab took two bytes @@ -554,33 +551,32 @@ function parse_line!( mapq % UInt8 end - # Fill in fields of the struct. The aux data begins at `i`. + # The aux data begins at `i`. # The integers have been checked to fit in an UInt32 by `parse_int`. - record.data = mem - record.qname_len = length(qname) % Int32 - record.tname_start = first(tname) % Int32 - record.tname_len = length(tname) % Int32 - record.aux_start = i % Int32 - record.qlen = qlen % UInt32 - record.qstart = qstart % UInt32 - record.qend = qend % UInt32 - record.tlen = tlen % UInt32 - record.tstart = tstart % UInt32 - record.tend = tend % UInt32 - record.matches = matches % UInt32 - record.alnlen = alnlen % UInt32 - record.mapq = mapq - record.strand = strand - - return record + return PAFRecord( + mem, + length(qname) % Int32, + first(tname) % Int32, + length(tname) % Int32, + i % Int32, + qlen % UInt32, + qstart % UInt32, + qend % UInt32, + tlen % UInt32, + tstart % UInt32, + tend % UInt32, + matches % UInt32, + alnlen % UInt32, + mapq, + strand, + ) end # Parse the remaining fields of an unmapped record. # These are validated like for mapped records, except that the integers may be zero, # as written by e.g. minimap2's `--paf-no-hit`. Their values are not stored. # i is the byte index of the strand plus two -function finish_unmapped!( - record::PAFRecord, +function finish_unmapped( mem::ImmutableMemoryView{UInt8}, qname::UnitRange{Int}, qlen::Int, @@ -597,24 +593,8 @@ function finish_unmapped!( (mapq, i) = @? parse_int(mem, i, true, true) mapq > 255 && return ParserException(i % Int32, Errors.IntegerOverflow) - # Fill in fields. The aux data begins at `i`. - record.data = mem - record.qname_len = length(qname) % Int32 - record.tname_start = 1 - record.tname_len = 0 - record.aux_start = i % Int32 - record.qlen = qlen % UInt32 - record.qstart = 0 - record.qend = 0 - record.tlen = 0 - record.tstart = 0 - record.tend = 0 - record.matches = 0 - record.alnlen = 0 - record.mapq = 0xff - record.strand = 0x00 - - return record + # The aux data begins at `i`. + return PAFRecord(mem, length(qname) % Int32, 1, 0, i % Int32, qlen % UInt32, 0, 0, 0, 0, 0, 0, 0, 0xff, 0x00) end # Just get the index of the next \t @@ -679,10 +659,11 @@ For efficiency, the auxiliary fields are not validated until they are accessed. If `io` is an `IO`, it is wrapped in a `BufferIO.BufReader` with a buffer size given by the keyword `buf_size`. -If `copy` is false, the *same* record will be returned on each iteration. -This record is a view into the reader's buffer, and is only valid until the reader is -advanced. This removes a few allocations and copies per iteration, but iterating the -reader will overwrite old records. +If `copy` is false, the iterated records are views into the reader's buffer, and are +only valid until the reader is advanced. Calling [`isempty`](@ref Base.isempty(::PAFReader)) +on the reader may also invalidate them. +This removes a few allocations and copies per iteration, but records must not be used +after the reader has moved on, e.g. by storing them. # Examples ```jldoctest @@ -702,30 +683,29 @@ julia> PAFReader(open(path_to_paf)) do reader 288659 julia> PAFReader(open(path_to_paf); copy=false) do reader - fst = first(reader) - all(reader) do record - record === fst - end + sum(record -> record.qlen, reader) end -true +889088 ``` See also: [`PAFRecord`](@ref), [`try_next!`](@ref) """ mutable struct PAFReader{I <: AbstractBufReader} const io::ExpandReader{I} - const record::PAFRecord # Number of bytes of the last parsed line, including its newline, which have not # yet been consumed from `io`. Consuming is deferred to the next call to `try_next!`, # such that the record, which is a view into the buffer of `io`, stays valid. to_consume::Int line::Int - # Whether to iterate copies of the record stored in this struct + # Whether to iterate copies of the parsed records copy::Bool + # Whether a record which is a view into the buffer of `io` has been returned to the + # user, and is still valid. While this is set, `isdone` must not fill the buffer. + view_live::Bool end function PAFReader(io::AbstractBufReader; copy::Bool = true) - return PAFReader{typeof(io)}(ExpandReader(io), PAFRecord(), 0, 1, copy) + return PAFReader{typeof(io)}(ExpandReader(io), 0, 1, copy, false) end function PAFReader(io::IO; buf_size::Int = 2^16, copy::Bool = true) @@ -766,12 +746,14 @@ function Base.copy(record::PAFRecord) end function Base.iterate(reader::PAFReader, (::Nothing) = nothing) - res = try_next!(reader) + res = next_view!(reader) return if res isa ParserException throw(res) elseif res === nothing nothing else + # A copy owns its data, so the view never reaches the user + reader.view_live = !reader.copy rec = reader.copy ? copy(res) : res (rec, nothing) end @@ -789,8 +771,9 @@ Try to read a record from the [`PAFReader`](@ref) and advance the reader. This m Even if the reader itself is not advanced, it may still fill its internal buffer, advancing the `IO` object it wraps. -The returned record is the reader's own record, which is a view into the reader's buffer. -It is overwritten, and its data invalidated, by the next call to `try_next!`. +The returned record is a view into the reader's buffer. Its data is invalidated by +the next call to `try_next!`, and may be invalidated by calling `isempty` on the reader. +Use `copy` on the record to obtain a record that owns its data. # Examples ```jldoctest @@ -815,12 +798,19 @@ true ``` """ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} + res = next_view!(reader) + res isa PAFRecord && (reader.view_live = true) + return res +end + +# Like `try_next!`, but the returned record is not yet exposed to the user +function next_view!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} consume_previous!(reader) (buffer, last_index, line_len) = buffer_line(reader.io) # If the line is empty, this is the end of the stream, if the rest of it is whitespace. iszero(last_index) && return check_tail(reader, buffer) # If this results in an error, we return the error. - res = parse_line!(reader.record, buffer[1:last_index]) + res = parse_line(buffer[1:last_index]) res isa ParserException && return ParserException(res.kind, getfield(res, :index_in_line), reader.line) # Else, we update the reader. This means `next!` will continuously return @@ -836,6 +826,7 @@ function consume_previous!(reader::PAFReader) consume(reader.io, reader.to_consume) reader.to_consume = 0 end + reader.view_live = false return nothing end @@ -884,29 +875,62 @@ function check_tail(reader::PAFReader, buffer::ImmutableMemoryView{UInt8})::Unio return end -# This is also used by `isempty`. It does not advance the reader, and does not -# invalidate any records returned previously. -function Base.isdone(reader::PAFReader, _...) +# This does not advance the reader, and never invalidates a returned record. +# If that is not possible, because a returned record is a view into the buffer and the +# buffer must be filled to determine the answer, return `missing`. +function Base.isdone(reader::PAFReader, _...)::Union{Bool, Missing} + buffer = get_buffer(reader.io) + rest = buffer[(reader.to_consume + 1):end] # Fast path: If the next line begins with any byte other than a newline, it is not # empty, so `try_next!` would return a record or an error. - buffer = get_buffer(reader.io) - n = reader.to_consume - if length(buffer) > n - b = @inbounds buffer[n + 1] + if !isempty(rest) + b = @inbounds rest[1] b ∈ (UInt8('\r'), UInt8('\n')) || return false end - # Slow path: We need to consume the previous line, and fill the buffer, which may - # overwrite the bytes the reader's record is a view into. So, make the record own its - # data first. This keeps any references to the record held by users valid. - if !iszero(n) - record = reader.record - setfield!(record, :data, copy(getfield(record, :data))) - consume_previous!(reader) + # If no view is live, we can freely consume the previous line and fill the buffer. + reader.view_live || return isdone_filling!(reader) + # Check if the next line is empty, mirroring `buffer_line`. + isempty(rest) && return missing + b = @inbounds rest[1] + if b == UInt8('\r') + length(rest) == 1 && return missing + @inbounds(rest[2]) == UInt8('\n') || return false + elseif b != UInt8('\n') + # Line is not empty, so `try_next!` returns a record or an error. + return false end + # The next line is empty, so the reader is done if the rest of the stream is + # whitespace, mirroring `check_tail`. We can only tell if it isn't. + for b in rest + b ∉ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) && return false + end + return missing +end + +function isdone_filling!(reader::PAFReader)::Bool + consume_previous!(reader) (buffer, last_index, _) = buffer_line(reader.io) return iszero(last_index) && check_tail(reader, buffer) === nothing end +""" + isempty(reader::PAFReader)::Bool + +Check whether `reader` has no more records, i.e. whether [`try_next!`](@ref) would return +`nothing`. This does not advance the reader, but may invalidate the previous record +by mutating the reader's inner state. + +In contrast, `Base.isdone(reader)`, never invalidates a record, +and returns `missing` if it can't determine whether the reader is done without invalidating +previous records. +""" +function Base.isempty(reader::PAFReader) + done = Base.isdone(reader) + done === missing || return done + # Invalidate the live view, such that `isdone` may fill the buffer + return isdone_filling!(reader) +end + @setup_workload begin path = joinpath(pkgdir(PairwiseMappingFormat), "test", "example.paf") @compile_workload begin diff --git a/test/runtests.jl b/test/runtests.jl index d6fab10..545dd2c 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -285,13 +285,17 @@ end records = PAFReader(collect, CursorReader(data)) @test [r.qname for r in records] == ["my_qname", "other", "my_qname"] - # With copy=false, the record is valid until the reader is advanced + # With copy=false, each record is a distinct view into the reader's buffer. + # CursorReader never moves its data, so the views stay valid here. PAFReader(CursorReader(data); copy = false) do reader rec = try_next!(reader) @test rec.qname == "my_qname" - @test try_next!(reader) === rec - @test rec.qname == "other" + rec2 = try_next!(reader) + @test rec2 !== rec + @test rec2.qname == "other" end + records = PAFReader(collect, CursorReader(data); copy = false) + @test [r.qname for r in records] == ["my_qname", "other", "my_qname"] # Parsing a vector aliases it, but copying the record does not v = collect(codeunits(good)) @@ -362,38 +366,72 @@ end @test length(collect(reader)) == 3 @test isempty(reader) - # isempty does not invalidate the record with copy=false, even when it needs to + # isempty does not invalidate records iterated with copy=true, even when it needs to # fill the buffer + expected = parse(PAFRecord, good) for chunks in ([good * "\n", good * "\n"], [good * "\n", "\n", "\r\n"], [good]) - reader = PAFReader(ChunkReader(chunks); copy = false) - rec = try_next!(reader) - expected = copy(rec) + reader = PAFReader(ChunkReader(chunks)) + (rec, _) = iterate(reader) while !isempty(reader) @test rec.qname == expected.qname @test rec.tname == expected.tname @test aux_data(rec) == aux_data(expected) - rec = try_next!(reader) - @test rec isa PAFRecord + (rec, _) = iterate(reader) end @test rec.qname == expected.qname - @test try_next!(reader) === nothing + @test rec.tname == expected.tname + @test iterate(reader) === nothing end - # isempty agrees with try_next! for all inputs and chunkings + # isdone does not fill the buffer when a view returned by try_next! is live, + # but isempty does + reader = PAFReader(ChunkReader([good * "\n", good * "\n"])) + rec = try_next!(reader) + @test Base.isdone(reader) === missing + @test rec.qname == "my_qname" + @test !isempty(reader) + @test try_next!(reader) isa PAFRecord + @test isempty(reader) + inputs = [ "", "\n", "\r\n", "\n\n \t\r\n", good, good * "\n", good * "\r\n\n", good * "\n\n" * good, good * "\nbad", good * "\n \n", good * "\n\r", - " ", "\nV", good * "\n\n\n\nx", + good * "\n\r\r\n", " ", "\nV", good * "\n\n\n\nx", ] - for data in inputs, split_at in 0:length(data) - # Empty chunks are not allowed, since fill_buffer returning zero signals EOF - chunks = filter!(!isempty, [data[1:split_at], data[(split_at + 1):end]]) + # Empty chunks are not allowed, since fill_buffer returning zero signals EOF + function chunkings(data) + splits = [filter!(!isempty, [data[1:i], data[(i + 1):end]]) for i in 0:length(data)] + return push!(splits, [string(c) for c in data]) + end + for data in inputs, chunks in chunkings(data) + # isdone and isempty agree with try_next!, and isdone never invalidates + # the record returned by try_next! reader = PAFReader(ChunkReader(chunks)) + (rec, expected) = (nothing, nothing) while true + maybe_done = Base.isdone(reader) + rec === nothing || @test getfield(rec, :data) == getfield(expected, :data) done = isempty(reader) + maybe_done === missing || @test maybe_done == done res = try_next!(reader) @test done == (res === nothing) res isa PAFRecord || break + (rec, expected) = (res, copy(res)) + end + + # When iterating with copy=true, no view is live, so isdone is never missing + reader = PAFReader(ChunkReader(chunks)) + while true + done = Base.isdone(reader) + @test done isa Bool + it = try + iterate(reader) + catch e + e isa PairwiseMappingFormat.ParserException || rethrow() + e + end + @test done == (it === nothing) + it isa Tuple || break end end end From 9f315cef0cbd2db62524159acaeb5441de0d6578 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 12:34:59 +0200 Subject: [PATCH 11/20] Have `try_next` respect reader's `copy` --- CHANGELOG.md | 5 ++++- docs/src/index.md | 11 +++++++---- src/PairwiseMappingFormat.jl | 34 +++++++++++++++++++--------------- test/runtests.jl | 20 +++++++++++++++----- 4 files changed, 45 insertions(+), 25 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 4cb322f..384c025 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -20,9 +20,12 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 of owning a copy of the data: - With `copy=false`, the record iterated by `PAFReader` and returned by `try_next!` is a view into the reader's buffer, and is invalidated when the reader is advanced. + Using an invalidated record is not detected. - `try_parse` and `parse` return a view into the input, so mutating the input also mutates the record. - `copy(record)` returns a record that owns its data. +- **Breaking:** `try_next!` now respects the `copy` keyword of `PAFReader`. With the default + `copy=true`, it returns a new record which owns its data. - **Breaking:** `PAFRecord` is now immutable. With `copy=false`, iterating a `PAFReader` no longer returns the same object on each iteration, but a new record each time. On Julia versions before 1.14, returning a `PAFRecord` from a non-inlined function, @@ -47,6 +50,6 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Fixed - `isempty(reader)` no longer consumes a record from a `PAFReader`. It may invalidate the most recently returned record if that record is a view into the reader's buffer, - i.e. if it was returned by `try_next!`, or iterated with `copy=false`. + i.e. if the reader was constructed with `copy=false`. `PAFReader` now implements `Base.isdone`, which never invalidates records, and returns `missing` when it can't determine whether the reader is done without doing so. diff --git a/docs/src/index.md b/docs/src/index.md index 70fcb4c..4158b44 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -47,13 +47,16 @@ julia> PAFReader(open(path_to_paf)) do reader ``` The [`PAFReader`](@ref) constructor takes an optional keyword `copy`, which defaults to `true`. -If it is `false`, the iterated [`PAFRecord`](@ref)s are views into the reader's buffer, -instead of owning a copy of their data. +If it is `false`, the [`PAFRecord`](@ref)s returned by iterating the reader, or by +[`try_next!`](@ref), are views into the reader's buffer, instead of owning a copy of their data. Such a record is only valid until the reader is advanced, and may be invalidated by calling [`isempty`](@ref Base.isempty(::PAFReader)) on the reader. This reduces allocations and copying and gives a slight speed increase, but may result -in bugs if records from previous iterations are stored, since their data may be -silently overwritten. Use `copy` on a record to obtain a record that owns its data: +in bugs if records from previous iterations are stored. +Using an invalidated record is not detected, and does not throw an error. Instead, its data +may be silently overwritten with the content of later lines, so it may appear to be a +different, or even malformed, record. +Use `copy` on a record to obtain a record that owns its data: ```jldoctest julia> PAFReader(open(path_to_paf); copy=false) do reader diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 5206548..b63ee0e 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -659,11 +659,13 @@ For efficiency, the auxiliary fields are not validated until they are accessed. If `io` is an `IO`, it is wrapped in a `BufferIO.BufReader` with a buffer size given by the keyword `buf_size`. -If `copy` is false, the iterated records are views into the reader's buffer, and are -only valid until the reader is advanced. Calling [`isempty`](@ref Base.isempty(::PAFReader)) -on the reader may also invalidate them. +If `copy` is false, the records returned by iterating the reader, or by [`try_next!`](@ref), +are views into the reader's buffer, and are only valid until the reader is advanced. +Calling [`isempty`](@ref Base.isempty(::PAFReader)) on the reader may also invalidate them. This removes a few allocations and copies per iteration, but records must not be used after the reader has moved on, e.g. by storing them. +Using an invalidated record is not detected: Its data may be silently overwritten +with the content of later lines. Use `copy` on a record to obtain a record that owns its data. # Examples ```jldoctest @@ -701,6 +703,7 @@ mutable struct PAFReader{I <: AbstractBufReader} copy::Bool # Whether a record which is a view into the buffer of `io` has been returned to the # user, and is still valid. While this is set, `isdone` must not fill the buffer. + # Only ever set when `copy` is false. view_live::Bool end @@ -746,16 +749,13 @@ function Base.copy(record::PAFRecord) end function Base.iterate(reader::PAFReader, (::Nothing) = nothing) - res = next_view!(reader) + res = try_next!(reader) return if res isa ParserException throw(res) elseif res === nothing nothing else - # A copy owns its data, so the view never reaches the user - reader.view_live = !reader.copy - rec = reader.copy ? copy(res) : res - (rec, nothing) + (res, nothing) end end @@ -771,9 +771,10 @@ Try to read a record from the [`PAFReader`](@ref) and advance the reader. This m Even if the reader itself is not advanced, it may still fill its internal buffer, advancing the `IO` object it wraps. -The returned record is a view into the reader's buffer. Its data is invalidated by -the next call to `try_next!`, and may be invalidated by calling `isempty` on the reader. -Use `copy` on the record to obtain a record that owns its data. +If the reader was constructed with `copy=true` (the default), the returned record owns +its data. If `copy=false`, the record is a view into the reader's buffer, which is +invalidated when the reader is advanced, and may be invalidated by calling +[`isempty`](@ref Base.isempty(::PAFReader)) on the reader. See [`PAFReader`](@ref). # Examples ```jldoctest @@ -799,11 +800,14 @@ true """ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} res = next_view!(reader) - res isa PAFRecord && (reader.view_live = true) + res isa PAFRecord || return res + # A copy owns its data, so the view never reaches the user + reader.copy && return copy(res) + reader.view_live = true return res end -# Like `try_next!`, but the returned record is not yet exposed to the user +# Like `try_next!`, but always returns a view, which is not yet exposed to the user function next_view!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} consume_previous!(reader) (buffer, last_index, line_len) = buffer_line(reader.io) @@ -917,8 +921,8 @@ end isempty(reader::PAFReader)::Bool Check whether `reader` has no more records, i.e. whether [`try_next!`](@ref) would return -`nothing`. This does not advance the reader, but may invalidate the previous record -by mutating the reader's inner state. +`nothing`. This does not advance the reader, but if the reader was constructed with +`copy=false`, it may invalidate the previous record by mutating the reader's inner state. In contrast, `Base.isdone(reader)`, never invalidates a record, and returns `missing` if it can't determine whether the reader is done without invalidating diff --git a/test/runtests.jl b/test/runtests.jl index 545dd2c..47cbfe4 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -297,6 +297,13 @@ end records = PAFReader(collect, CursorReader(data); copy = false) @test [r.qname for r in records] == ["my_qname", "other", "my_qname"] + # try_next! respects copy: With copy=true, records stay valid after advancing, + # even when the reader's buffer is overwritten + reader = PAFReader(ChunkReader([good * "\n", other * "\n"])) + rec = try_next!(reader) + @test try_next!(reader).qname == "other" + @test rec.qname == "my_qname" + # Parsing a vector aliases it, but copying the record does not v = collect(codeunits(good)) rec = try_parse(v) @@ -385,7 +392,7 @@ end # isdone does not fill the buffer when a view returned by try_next! is live, # but isempty does - reader = PAFReader(ChunkReader([good * "\n", good * "\n"])) + reader = PAFReader(ChunkReader([good * "\n", good * "\n"]); copy = false) rec = try_next!(reader) @test Base.isdone(reader) === missing @test rec.qname == "my_qname" @@ -403,20 +410,23 @@ end splits = [filter!(!isempty, [data[1:i], data[(i + 1):end]]) for i in 0:length(data)] return push!(splits, [string(c) for c in data]) end - for data in inputs, chunks in chunkings(data) + for data in inputs, chunks in chunkings(data), copy in (true, false) # isdone and isempty agree with try_next!, and isdone never invalidates - # the record returned by try_next! - reader = PAFReader(ChunkReader(chunks)) + # the record returned by try_next!. With copy=true, no view is live, so + # isdone is never missing, and isempty does not invalidate the record either. + reader = PAFReader(ChunkReader(chunks); copy) (rec, expected) = (nothing, nothing) while true maybe_done = Base.isdone(reader) + copy && @test maybe_done isa Bool rec === nothing || @test getfield(rec, :data) == getfield(expected, :data) done = isempty(reader) + copy && rec !== nothing && @test getfield(rec, :data) == getfield(expected, :data) maybe_done === missing || @test maybe_done == done res = try_next!(reader) @test done == (res === nothing) res isa PAFRecord || break - (rec, expected) = (res, copy(res)) + (rec, expected) = (res, Base.copy(res)) end # When iterating with copy=true, no view is live, so isdone is never missing From c0761dfb46b2697797319abc448cef501f9fb1d3 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 13:27:37 +0200 Subject: [PATCH 12/20] Apply whitespace rules more consistently --- CHANGELOG.md | 12 ++++ docs/src/index.md | 12 ++-- src/PairwiseMappingFormat.jl | 115 +++++++++++++++++++++++------------ test/runtests.jl | 78 ++++++++++++++++++++++-- 4 files changed, 167 insertions(+), 50 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 384c025..7a1917c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -44,6 +44,18 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - **Breaking:** Fields 6-12 of unmapped records are now validated like those of mapped records, except that the integer fields may be zero. Previously, their content was ignored. +- `parse(PAFRecord, x)` now accepts the same inputs as `try_parse`, i.e. any `x` that + implements `MemoryView`, such as `Vector{UInt8}`. Previously, it only accepted + `String` and `SubString{String}`. +- **Breaking:** `parse` and `try_parse` now strip one optional trailing `\n` or `\r\n`, + and reject any other `\n` with the new error kind `EmbeddedNewline`. + Previously, a trailing newline was part of the last field, and an embedded newline + was accepted in names and auxiliary data. +- **Breaking:** In `PAFReader`, the first line consisting only of whitespace now marks the + end of the records, after which all input must be whitespace. + Previously, only an empty line did, and a line of e.g. a single space was a + `TooFewFields` error. + ### Removed - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. diff --git a/docs/src/index.md b/docs/src/index.md index 4158b44..12c3a43 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -273,11 +273,12 @@ true ``` ## Misc info -* The `PAFReader` cannot handle trailing whitespace, except trailing newlines at - the end of the file. This is because trailing whitespace may be significant - in some records, e.g. if it ends with an auxiliary field of the element type - `Z`. - A PAF record with trailing whitespace is considered invalid. +* Lines end with `\n` or `\r\n`. Other whitespace at the end of a record is not + stripped, because it may be significant, e.g. if the record ends with an auxiliary + field of the element type `Z`. Instead, it is part of the last field, and may make it + invalid. +* The first line in a file consisting only of whitespace marks the end of the records. + All remaining input must then be whitespace. ## Reference @@ -294,6 +295,7 @@ Base.isempty(::PAFReader) is_mapped PairwiseMappingFormat.Errors PairwiseMappingFormat.try_parse +Base.parse(::Type{PAFRecord}, ::Any) PairwiseMappingFormat.ParserException PairwiseMappingFormat.UnmappedError ``` diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index b63ee0e..37507ce 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -380,6 +380,7 @@ module Errors InvalidZero PositionOutOfBounds LineTooLong + EmbeddedNewline end end # module Errors @@ -451,6 +452,8 @@ function Base.showerror(io::IO, err::ParserException) "Query or target start or end position is longer than query/target length" elseif kind == Errors.LineTooLong "Line is too long. PairwiseMappingFormat.jl can't handle lines of 2147483647 bytes or longer" + elseif kind == Errors.EmbeddedNewline + "Line contains a newline. Only a single trailing \"\\n\" or \"\\r\\n\" is allowed" end print(buf, s) return print(io, String(take!(buf))) @@ -471,6 +474,11 @@ end Try parsing `x` into a `PAFRecord`. `x` may be any type that implements `MemoryView`, such as `String` or `Vector{UInt8}`. +`x` must contain a single line, optionally ending with `"\\n"` or `"\\r\\n"`. +Any other `\\n` in `x` gives an error of kind `EmbeddedNewline`. +Other trailing whitespace, including a lone `\\r`, is not stripped: +It is part of the last field, and may make it invalid. + !!! warning The resulting record is a view into the memory of `x`, and does not copy it. If `x` is mutated, the record is mutated too, and may become invalid. @@ -487,13 +495,47 @@ julia> typeof(PAF.try_parse("invalid string")) PairwiseMappingFormat.ParserException ``` -See also: [`PAFRecord`](@ref), [`try_next!`](@ref) +See also: [`PAFRecord`](@ref), [`try_next!`](@ref), [`parse`](@ref Base.parse(::Type{PAFRecord}, ::Any)) """ try_parse(x) = try_parse(ImmutableMemoryView(x)) -try_parse(x::ImmutableMemoryView{UInt8}) = parse_line(x) -function Base.parse(::Type{PAFRecord}, s::Union{String, SubString{String}}) - y = try_parse(s) +# Strip one optional trailing \n or \r\n +function strip_newline(x::ImmutableMemoryView{UInt8})::ImmutableMemoryView{UInt8} + (isempty(x) || @inbounds(x[end]) != UInt8('\n')) && return x + has_cr = length(x) > 1 && @inbounds(x[end - 1]) == UInt8('\r') + return @inbounds x[1:(end - 1 - has_cr)] +end + +function try_parse(x::ImmutableMemoryView{UInt8}) + x = strip_newline(x) + newline = findfirst(==(UInt8('\n')), x) + newline === nothing || return ParserException(newline, Errors.EmbeddedNewline) + return parse_line(x) +end + +""" + parse(::Type{PAFRecord}, x)::PAFRecord + +Parse `x` into a `PAFRecord`, or throw a [`ParserException`](@ref +PairwiseMappingFormat.ParserException) if `x` is not a valid PAF line. +`x` may be any type that implements `MemoryView`, such as `String` or `Vector{UInt8}`. +See [`try_parse`](@ref PairwiseMappingFormat.try_parse) for what counts as a line. + +!!! warning + The resulting record is a view into the memory of `x`, and does not copy it. + If `x` is mutated, the record is mutated too, and may become invalid. + Use `copy` on the record to obtain a record that owns its data. + +# Examples +```jldoctest +julia> record = parse(PAFRecord, codeunits("q\\t10\\t0\\t5\\t+\\tt\\t20\\t5\\t10\\t5\\t5\\t60\\n")); + +julia> (record.qname, record.tname, record.mapq) +("q", "t", 60) +``` +""" +function Base.parse(::Type{PAFRecord}, x) + y = try_parse(x) return y isa PAFRecord ? y : throw(y) end @@ -656,6 +698,12 @@ Iterating this object returns a `PAFRecord` for each valid line of input, and throws a [`ParserException`](@ref) when an invalid record is found. For efficiency, the auxiliary fields are not validated until they are accessed. +Lines end with `"\\n"` or `"\\r\\n"`. Other whitespace at the end of a line is not stripped: +It is part of the last field, and may make it invalid. +The first line consisting only of whitespace (spaces, tabs, `\\r` and `\\n`) marks the end +of the records. All remaining input must then be whitespace, else a `ParserException` of +kind `EmptyNonTailingLine` is thrown. + If `io` is an `IO`, it is wrapped in a `BufferIO.BufReader` with a buffer size given by the keyword `buf_size`. @@ -811,12 +859,15 @@ end function next_view!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} consume_previous!(reader) (buffer, last_index, line_len) = buffer_line(reader.io) - # If the line is empty, this is the end of the stream, if the rest of it is whitespace. - iszero(last_index) && return check_tail(reader, buffer) - # If this results in an error, we return the error. - res = parse_line(buffer[1:last_index]) - res isa ParserException && + line = @inbounds buffer[1:last_index] + res = parse_line(line) + if res isa ParserException + # A whitespace-only line, which can never be a valid record, is the end of the + # stream if the rest of it is whitespace. Checking only on errors keeps the + # check off the hot path. + all(is_whitespace, line) && return check_tail(reader, buffer) return ParserException(res.kind, getfield(res, :index_in_line), reader.line) + end # Else, we update the reader. This means `next!` will continuously return # an error if an error is found. reader.to_consume = line_len @@ -855,9 +906,11 @@ function buffer_line(io::ExpandReader)::Tuple{ImmutableMemoryView{UInt8}, Int, I return (buffer, newline - (1 + has_windows_newline), newline) end -# Called when the next line is empty. Return `nothing` if the rest of the stream is -# whitespace, else an error. -# The buffer may end after the empty line, so keep filling it until EOF. +is_whitespace(b::UInt8) = b ∈ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) + +# Called when the next line is whitespace-only, with `buffer` beginning at that line. +# Return `nothing` if the rest of the stream is whitespace, else an error. +# The buffer may end after the line, so keep filling it until EOF. function check_tail(reader::PAFReader, buffer::ImmutableMemoryView{UInt8})::Union{Nothing, ParserException} io = reader.io line = reader.line @@ -868,7 +921,7 @@ function check_tail(reader::PAFReader, buffer::ImmutableMemoryView{UInt8})::Unio col += 1 line += b == UInt8('\n') col *= b != UInt8('\n') - if b ∉ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) + if !is_whitespace(b) return ParserException(Errors.EmptyNonTailingLine, col % Int32, line) end end @@ -883,38 +936,20 @@ end # If that is not possible, because a returned record is a view into the buffer and the # buffer must be filled to determine the answer, return `missing`. function Base.isdone(reader::PAFReader, _...)::Union{Bool, Missing} - buffer = get_buffer(reader.io) - rest = buffer[(reader.to_consume + 1):end] - # Fast path: If the next line begins with any byte other than a newline, it is not - # empty, so `try_next!` would return a record or an error. - if !isempty(rest) - b = @inbounds rest[1] - b ∈ (UInt8('\r'), UInt8('\n')) || return false - end + rest = get_buffer(reader.io)[(reader.to_consume + 1):end] + # The reader is done exactly when the rest of the stream is whitespace: + # If the next line contains non-whitespace, `try_next!` returns a record or an error. + # Else, it is the end of the stream only if the rest of it is whitespace. + # This usually returns at the first byte. + all(is_whitespace, rest) || return false # If no view is live, we can freely consume the previous line and fill the buffer. - reader.view_live || return isdone_filling!(reader) - # Check if the next line is empty, mirroring `buffer_line`. - isempty(rest) && return missing - b = @inbounds rest[1] - if b == UInt8('\r') - length(rest) == 1 && return missing - @inbounds(rest[2]) == UInt8('\n') || return false - elseif b != UInt8('\n') - # Line is not empty, so `try_next!` returns a record or an error. - return false - end - # The next line is empty, so the reader is done if the rest of the stream is - # whitespace, mirroring `check_tail`. We can only tell if it isn't. - for b in rest - b ∉ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) && return false - end - return missing + # Else, we can't tell if the stream continues with non-whitespace. + return reader.view_live ? missing : isdone_filling!(reader) end function isdone_filling!(reader::PAFReader)::Bool consume_previous!(reader) - (buffer, last_index, _) = buffer_line(reader.io) - return iszero(last_index) && check_tail(reader, buffer) === nothing + return check_tail(reader, get_buffer(reader.io)) === nothing end """ diff --git a/test/runtests.jl b/test/runtests.jl index 47cbfe4..caf28ca 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -329,11 +329,48 @@ end test_aux_content(good * stringaux_with_whitespace, goodaux * stringaux_with_whitespace) end -# Note: We currently DONT trim whitespace off ends of records, for the sole reason that -# whitespace is allowed in string AuxData. -# This is very annoying, because it makes parsing less robust to erroneous trailing whitespace. -# However, I can't think of a good way to deal with this. -# TODO: Find a better method +@testset "Parsing lines" begin + good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" + expected = parse(PAFRecord, good) + same(rec) = rec isa PAFRecord && getfield(rec, :data) == getfield(expected, :data) + + # parse accepts the same inputs as try_parse + for x in (good, SubString(good, 1), codeunits(good), collect(codeunits(good)), MemoryView(good)) + @test same(parse(PAFRecord, x)) + @test same(try_parse(x)) + end + @test_throws PairwiseMappingFormat.ParserException parse(PAFRecord, codeunits("abc")) + + # One optional trailing newline is stripped, also before aux data + @test same(try_parse(good * "\n")) + @test same(try_parse(good * "\r\n")) + @test aux_data(parse(PAFRecord, good * "\tab:Z:x\r\n"))["ab"] == "x" + + # Any other newline is an error + @test try_parse(good * "\n\n").kind == Errors.EmbeddedNewline + @test try_parse(good * "\r\n\n").kind == Errors.EmbeddedNewline + @test try_parse(good * "\n" * good).kind == Errors.EmbeddedNewline + @test try_parse(replace(good, "my_qname" => "my\nqname")).kind == Errors.EmbeddedNewline + @test try_parse(good * "\tab:Z:x\ny").kind == Errors.EmbeddedNewline + @test try_parse("\n").kind == Errors.TooFewFields + + # An embedded \r is data + @test parse(PAFRecord, replace(good, "my_qname" => "my\rqname")).qname == "my\rqname" + + # Other trailing whitespace is not stripped, but part of the last field + @test try_parse(good * " ").kind == Errors.BadInteger + @test try_parse(good * "\r").kind == Errors.BadInteger + # A trailing tab begins an empty aux field, which is only invalid when accessed + @test try_parse(good * "\t") isa PAFRecord + @test try_parse(good * " \n").kind == Errors.BadInteger + @test aux_data(parse(PAFRecord, good * "\tab:Z:x \r\n"))["ab"] == "x " + # A trailing \r is kept in aux data, where it is invalid for Z values + @test collect(MemoryView(aux_data(parse(PAFRecord, good * "\tab:Z:x\r")))) == codeunits("ab:Z:x\r") +end + +# Note: We don't trim whitespace off ends of records, because whitespace is allowed in +# string AuxData. Trailing whitespace is part of the last field, and may make it invalid. +# After the records, a file may end with whitespace-only lines. @testset "Trailing whitespace" begin data = "my_qname 301156 3 301145 + my_tname 6701780 2764860 3066002 301142 301142 0 an:Z: dljdl \r\n" data *= "some 5 3 4 + other 9 5 8 2 2 0\r\n" @@ -361,6 +398,36 @@ end # Trailing whitespace delivered in multiple chunks is still fine reader = PAFReader(ChunkReader([good * "\n\n", "\r\n", "\n"])) @test length(collect(reader)) == 1 + + # The first whitespace-only line ends the records, whether or not it is empty + for tail in ("\n \n", "\n \n\n", "\n\t\r\n \n", "\n ", "\n\r", "\n\r\r\n", "\n\n \n") + @test length(PAFReader(collect, IOBuffer(good * tail))) == 1 + reader = PAFReader(IOBuffer(good * tail)) + @test try_next!(reader) isa PAFRecord + @test isempty(reader) + @test try_next!(reader) === nothing + end + @test isempty(PAFReader(IOBuffer(" \n\t\n"))) + + # ... after which only whitespace may follow + for tail in ("\n \n" * good, "\n\t\nx", "\n \r\n\n" * good * "\n") + reader = PAFReader(IOBuffer(good * tail)) + @test try_next!(reader) isa PAFRecord + @test !isempty(reader) + @test try_next!(reader).kind == Errors.EmptyNonTailingLine + end + + # Trailing whitespace on a record is not stripped + reader = PAFReader(IOBuffer(good * " \n")) + @test try_next!(reader).kind == Errors.BadInteger + reader = PAFReader(IOBuffer(good * "\tab:Z:x \n" * good * "\t\n")) + @test aux_data(try_next!(reader))["ab"] == "x " + @test try_next!(reader) isa PAFRecord # empty aux field is only invalid when accessed + + # A line with non-whitespace after leading whitespace is an ordinary error + reader = PAFReader(IOBuffer(good * "\n \tx\n")) + @test try_next!(reader) isa PAFRecord + @test try_next!(reader).kind == Errors.BadInteger end @testset "isempty" begin @@ -404,6 +471,7 @@ end "", "\n", "\r\n", "\n\n \t\r\n", good, good * "\n", good * "\r\n\n", good * "\n\n" * good, good * "\nbad", good * "\n \n", good * "\n\r", good * "\n\r\r\n", " ", "\nV", good * "\n\n\n\nx", + good * "\n \n" * good, good * "\n\t \r\n ", " \tx", good * "\n \n\tx", ] # Empty chunks are not allowed, since fill_buffer returning zero signals EOF function chunkings(data) From e51cbae31a421dcc8f614165b9592a2e12683700 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 13:40:37 +0200 Subject: [PATCH 13/20] Make malformed records advance the reader also --- CHANGELOG.md | 7 +++ docs/src/index.md | 21 ++++---- src/PairwiseMappingFormat.jl | 94 +++++++++++++++++++++--------------- test/runtests.jl | 70 +++++++++++++++++++++++++-- 4 files changed, 138 insertions(+), 54 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 7a1917c..459f5eb 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -56,6 +56,13 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 Previously, only an empty line did, and a line of e.g. a single space was a `TooFewFields` error. +- **Breaking:** When `try_next!` returns a `ParserException`, or iterating a `PAFReader` + throws one, the reader now advances past the invalid line, so the next call continues + from the following line. Previously, the reader returned the same error forever. + This allows skipping invalid lines. For an `EmptyNonTailingLine` error, the reader + skips all consecutive whitespace-only lines, and the error is reported at the first + of them, instead of at the first non-whitespace byte after them. + ### Removed - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. diff --git a/docs/src/index.md b/docs/src/index.md index 12c3a43..a5cf17f 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -189,13 +189,11 @@ PairwiseMappingFormat.ParserException(PairwiseMappingFormat.Errors.TooFewFields, Similarly, the next record of a `PAFReader` may be obtained with the unexported [`try_next!`](@ref) function. This function will either: -* Return `nothing` if the reader has no more records -* Return a `ParserException` if the next record is invalid -* Return a `PAFRecord` if there is another valid record. - Only in this case it will advance the position of the reader. - -In other words, if a call to [`try_next!`](@ref) returns `nothing` or a `ParserException`, -any future calls will return the same value forever: +* Return `nothing` if the reader has no more records. Any future calls will + also return `nothing`. +* Return a `ParserException` if the next record is invalid, and advance the reader + past the invalid line. Hence, invalid lines can be skipped by ignoring the error. +* Return a `PAFRecord` if there is another valid record, and advance the reader. ```jldoctest public; output = false reader = PAFReader(open(path_to_paf)) @@ -210,12 +208,13 @@ nothings = [try_next!(reader) for i in 1:10] close(reader) -reader = PAFReader(IOBuffer("not a PAF record")) +reader = PAFReader(IOBuffer("not a PAF record\n" * open(readline, path_to_paf))) err = try_next!(reader) - @assert err isa PAF.ParserException -err2 = try_next!(reader) -@assert err == err2 + +# The invalid line was skipped, so the next record can be read +@assert try_next!(reader) isa PAFRecord +@assert try_next!(reader) === nothing # output diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 37507ce..7e24b10 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -696,24 +696,14 @@ end Construct a `PAFReader`, an iterator of [`PAFRecord`](@ref) read from `io`. Iterating this object returns a `PAFRecord` for each valid line of input, and throws a [`ParserException`](@ref) when an invalid record is found. -For efficiency, the auxiliary fields are not validated until they are accessed. - -Lines end with `"\\n"` or `"\\r\\n"`. Other whitespace at the end of a line is not stripped: -It is part of the last field, and may make it invalid. -The first line consisting only of whitespace (spaces, tabs, `\\r` and `\\n`) marks the end -of the records. All remaining input must then be whitespace, else a `ParserException` of -kind `EmptyNonTailingLine` is thrown. - -If `io` is an `IO`, it is wrapped in a `BufferIO.BufReader` with a buffer size -given by the keyword `buf_size`. If `copy` is false, the records returned by iterating the reader, or by [`try_next!`](@ref), -are views into the reader's buffer, and are only valid until the reader is advanced. -Calling [`isempty`](@ref Base.isempty(::PAFReader)) on the reader may also invalidate them. +are views into the reader's buffer, and are only valid until the reader is mutated, +which happens on iterating the reader or, `try_next!`, `isempty`. This removes a few allocations and copies per iteration, but records must not be used -after the reader has moved on, e.g. by storing them. -Using an invalidated record is not detected: Its data may be silently overwritten -with the content of later lines. Use `copy` on a record to obtain a record that owns its data. +after they are invalidated. + +See also: [`PAFRecord`](@ref), [`try_next!`](@ref) # Examples ```jldoctest @@ -738,7 +728,15 @@ julia> PAFReader(open(path_to_paf); copy=false) do reader 889088 ``` -See also: [`PAFRecord`](@ref), [`try_next!`](@ref) +# Extended help +After throwing a `ParserException`, the reader advances over the malformed record, +so iteration may be consumed after catching the exception. + +For efficiency, the auxiliary fields are not validated until they are accessed. + +Lines end with `"\\n"` or `"\\r\\n"`. Other whitespace at the end of a line is not stripped: +It is part of the last field, and may make it invalid. + """ mutable struct PAFReader{I <: AbstractBufReader} const io::ExpandReader{I} @@ -813,12 +811,10 @@ end Try to read a record from the [`PAFReader`](@ref) and advance the reader. This may yield one of three options: * If the reader is empty, return `nothing` and do not advance the reader * If the next record is invalid, return (do not throw) a `ParserException` and - do not advance the reader + advance the reader past the invalid line. The next call continues from the next line, + so invalid lines may be skipped by ignoring the error. * If the next record is valid, return it as a [`PAFRecord`](@ref) and advance the reader -Even if the reader itself is not advanced, it may still fill its internal buffer, advancing -the `IO` object it wraps. - If the reader was constructed with `copy=true` (the default), the returned record owns its data. If `copy=false`, the record is a view into the reader's buffer, which is invalidated when the reader is advanced, and may be invalidated by calling @@ -837,12 +833,15 @@ true julia> typeof([try_next!(reader) for i in 1:100]) Vector{Nothing} (alias for Array{Nothing, 1}) -julia> close(reader); reader = PAFReader(IOBuffer("bad data")); +julia> close(reader); reader = PAFReader(IOBuffer("bad data\\n" * open(readline, path_to_paf))); julia> try_next!(reader) isa PAF.ParserException true -julia> all(i -> try_next!(reader) isa PAF.ParserException, 1:100) +julia> try_next!(reader) isa PAFRecord +true + +julia> try_next!(reader) === nothing true ``` """ @@ -861,15 +860,26 @@ function next_view!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserExceptio (buffer, last_index, line_len) = buffer_line(reader.io) line = @inbounds buffer[1:last_index] res = parse_line(line) - if res isa ParserException - # A whitespace-only line, which can never be a valid record, is the end of the - # stream if the rest of it is whitespace. Checking only on errors keeps the - # check off the hot path. - all(is_whitespace, line) && return check_tail(reader, buffer) - return ParserException(res.kind, getfield(res, :index_in_line), reader.line) + # A whitespace-only line, which can never be a valid record, is the end of the + # stream if the rest of it is whitespace. Checking only on errors keeps the + # check off the hot path. + if res isa ParserException && all(is_whitespace, line) + tail = check_tail(reader, buffer) + tail === nothing && return nothing + # Skip all the whitespace-only lines, and report the error at the first of them + (n_bytes, n_lines) = tail + err = ParserException(Errors.EmptyNonTailingLine, Int32(1), reader.line) + reader.to_consume = n_bytes + reader.line += n_lines + return err + end + res = if res isa ParserException + ParserException(res.kind, getfield(res, :index_in_line), reader.line) + else + res end - # Else, we update the reader. This means `next!` will continuously return - # an error if an error is found. + # Advance past the line, even if it is invalid, such that the next call + # continues from the next line reader.to_consume = line_len reader.line += 1 return res @@ -909,20 +919,24 @@ end is_whitespace(b::UInt8) = b ∈ (UInt8('\r'), UInt8('\n'), UInt8('\t'), UInt8(' ')) # Called when the next line is whitespace-only, with `buffer` beginning at that line. -# Return `nothing` if the rest of the stream is whitespace, else an error. +# Return `nothing` if the rest of the stream is whitespace. Else, return the number of +# bytes and lines of the whitespace-only lines before the first line with non-whitespace. # The buffer may end after the line, so keep filling it until EOF. -function check_tail(reader::PAFReader, buffer::ImmutableMemoryView{UInt8})::Union{Nothing, ParserException} +# This does not mutate `reader`, except by filling its buffer. +function check_tail(reader::PAFReader, buffer::ImmutableMemoryView{UInt8})::Union{Nothing, Tuple{Int, Int}} io = reader.io - line = reader.line - col = 0 + n_lines = 0 + # Index of the first byte after the last newline seen + line_start = 1 n_scanned = 0 while true - @inbounds for b in buffer[(n_scanned + 1):end] - col += 1 - line += b == UInt8('\n') - col *= b != UInt8('\n') - if !is_whitespace(b) - return ParserException(Errors.EmptyNonTailingLine, col % Int32, line) + for i in (n_scanned + 1):length(buffer) + b = @inbounds buffer[i] + if b == UInt8('\n') + n_lines += 1 + line_start = i + 1 + elseif !is_whitespace(b) + return (line_start - 1, n_lines) end end n_scanned = length(buffer) diff --git a/test/runtests.jl b/test/runtests.jl index caf28ca..c14b6ae 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -430,6 +430,51 @@ end @test try_next!(reader).kind == Errors.BadInteger end +@testset "Skipping invalid lines" begin + good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" + kinds(results) = [r isa PAFRecord ? :record : (r.kind, r.line) for r in results] + function read_all(reader) + results = [] + while (res = try_next!(reader)) !== nothing + push!(results, res) + end + return results + end + + # try_next! advances past invalid lines + data = join([good, "bad", good, good * "\tx", good * " ", good], '\n') + @test kinds(read_all(PAFReader(IOBuffer(data)))) == + [:record, (Errors.TooFewFields, 2), :record, :record, (Errors.BadInteger, 5), :record] + + # Including a truncated last line, as from a killed job + reader = PAFReader(IOBuffer(good * "\n" * good[1:30])) + @test try_next!(reader) isa PAFRecord + @test try_next!(reader).line == 2 + @test try_next!(reader) === nothing + @test isempty(reader) + + # A run of whitespace-only lines followed by non-whitespace is one error, + # reported at the first whitespace-only line. Reading continues after the run. + data = good * "\n\n \r\n\t\n" * good * "\n\nbad\n\n" + @test kinds(read_all(PAFReader(IOBuffer(data)))) == [ + :record, (Errors.EmptyNonTailingLine, 2), :record, + (Errors.EmptyNonTailingLine, 6), (Errors.TooFewFields, 7), + ] + + # Iteration throws, but also advances, so it can be resumed after catching + reader = PAFReader(IOBuffer(join([good, "bad", good], '\n'))) + @test iterate(reader) isa Tuple + @test_throws PairwiseMappingFormat.ParserException iterate(reader) + @test iterate(reader) isa Tuple + @test iterate(reader) === nothing + + # Mixing iteration and try_next! gives the same results + reader = PAFReader(IOBuffer(join(["bad", "bad", good], '\n'))) + @test_throws PairwiseMappingFormat.ParserException iterate(reader) + @test try_next!(reader).line == 2 + @test only(collect(reader)).qname == "my_qname" +end + @testset "isempty" begin good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" @@ -472,16 +517,35 @@ end good * "\n\n" * good, good * "\nbad", good * "\n \n", good * "\n\r", good * "\n\r\r\n", " ", "\nV", good * "\n\n\n\nx", good * "\n \n" * good, good * "\n\t \r\n ", " \tx", good * "\n \n\tx", + "bad\n" * good, good * "\nbad\n" * good * "\n", "bad\r\nbad\n\n \n" * good, ] # Empty chunks are not allowed, since fill_buffer returning zero signals EOF function chunkings(data) splits = [filter!(!isempty, [data[1:i], data[(i + 1):end]]) for i in 0:length(data)] return push!(splits, [string(c) for c in data]) end + # Summarize the results of reading all of `reader`, continuing past errors. + # Each call must advance the reader, so the number of results is bounded. + function read_all(reader, max_results) + results = [] + while (res = try_next!(reader)) !== nothing + push!(results, res isa PAFRecord ? copy(getfield(res, :data)) : res) + length(results) > max_results && error("Reader did not advance") + end + return results + end + for data in inputs + expected_results = read_all(PAFReader(IOBuffer(data)), ncodeunits(data)) + for chunks in chunkings(data) + @test read_all(PAFReader(ChunkReader(chunks)), ncodeunits(data)) == expected_results + end + end + for data in inputs, chunks in chunkings(data), copy in (true, false) # isdone and isempty agree with try_next!, and isdone never invalidates # the record returned by try_next!. With copy=true, no view is live, so # isdone is never missing, and isempty does not invalidate the record either. + # Reading continues past errors. reader = PAFReader(ChunkReader(chunks); copy) (rec, expected) = (nothing, nothing) while true @@ -493,8 +557,8 @@ end maybe_done === missing || @test maybe_done == done res = try_next!(reader) @test done == (res === nothing) - res isa PAFRecord || break - (rec, expected) = (res, Base.copy(res)) + res === nothing && break + (rec, expected) = res isa PAFRecord ? (res, Base.copy(res)) : (nothing, nothing) end # When iterating with copy=true, no view is live, so isdone is never missing @@ -509,7 +573,7 @@ end e end @test done == (it === nothing) - it isa Tuple || break + it === nothing && break end end end From 13da3a813264bcf5357e044200548592255583e6 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 13:49:56 +0200 Subject: [PATCH 14/20] Small documentation and typo fixes --- docs/src/index.md | 21 ++++++++++++++------- src/PairwiseMappingFormat.jl | 26 +++++++++++++++++++++----- 2 files changed, 35 insertions(+), 12 deletions(-) diff --git a/docs/src/index.md b/docs/src/index.md index a5cf17f..0a82629 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -105,7 +105,7 @@ underlying fields as they are stored in the `PAFRecord`. in the query and target, respectively. This does not take strandedness into account. * `qend` and `tend`. Of type `Int`. Rightmost (i.e. highest) position of the alignment in the query and target, respectively. This does not take strandedness into account. - As such, we always have `record.qend - record.qtart) >= 0`, and likewise for the target. + As such, we always have `record.qend - record.qstart >= 0`, and likewise for the target. * `alnlen::Int`. Length of alignment. That includes gaps in query and target, and may therefore not match either `qend - qstart + 1`, or that of the target. * `matches::Int`. Number of residue matches (not mismatches) in the alignment. @@ -118,7 +118,7 @@ underlying fields as they are stored in the `PAFRecord`. The PAF format does not contain fields for alignment identity, query coverage or target coverage. However, these can be approximated with the functions [`aln_identity`](@ref), -[`query_coverage`](@ref) and [`query_coverage`](@ref). +[`query_coverage`](@ref) and [`target_coverage`](@ref). Besides being iterated from [`PAFReader`](@ref), records can also be created by parsing from a bytes-like object such as a `String`: @@ -171,7 +171,7 @@ ERROR: Error when parsing PAF record on line 1, near byte number 14 in line: Not [...] ``` -In some circuomstances, throwing exceptions may not be acceptable, and so PairwiseMappingFormat.jl +In some circumstances, throwing exceptions may not be acceptable, and so PairwiseMappingFormat.jl provide functionality for returning errors as values. These values can then be checked to programmatically handle error states. @@ -186,7 +186,7 @@ julia> println(PAF.try_parse("not a PAF line")) PairwiseMappingFormat.ParserException(PairwiseMappingFormat.Errors.TooFewFields, 14, 1) ``` -Similarly, the next record of a `PAFReader` may be obtained with the unexported [`try_next!`](@ref) +Similarly, the next record of a `PAFReader` may be obtained with the [`try_next!`](@ref) function. This function will either: * Return `nothing` if the reader has no more records. Any future calls will @@ -231,9 +231,16 @@ julia> err.kind TooFewFields::Err = 0 ``` -The precise value of this `.kind` object for any given error condition is subject -to change across minor versions and new values may be introduced. -This is because the same error may be detected in multiple different ways. +The `.kind` is intended for diagnostics, such as logging or counting errors, and should +not be relied on to handle errors programmatically. +An invalid line is often invalid in several ways, and which of them is reported depends +on implementation details. Therefore, the kind reported for a given invalid input is +deterministic within a version of this package, but may change in any minor release, +and new kinds may be added. +In contrast, whether input is rejected, and the line the error is reported on, are stable. + +The error message also includes the approximate position of the error within the line. +This position is only meant for humans to read, and is not part of the API. ## [Auxiliary fields](@id aux) [`PAFRecord`](@ref)s may contain extra auxiliary fields at the end of the records, diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 7e24b10..81d59c5 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -1,7 +1,7 @@ """ module PairwiseMappingFormat -Parse files of the PairwisemAppingFormat (PAF) format. +Parse files in the Pairwise mApping Format (PAF). """ module PairwiseMappingFormat @@ -97,7 +97,7 @@ struct PAFRecord matches::UInt32 alnlen::UInt32 # always >= 1 mapq::UInt8 # stored as 0xff for missing - strand::UInt8 # 0x00: *, 0x01: -, 0x02: + + strand::UInt8 # 0x00: *, 0x01: +, 0x02: - # 16 bits of padding, which we can use for something else later end @@ -358,6 +358,10 @@ end Wrapper module used to contain the instances of error types, and not pollute the namespace of the package. +The error kinds are intended for diagnostics. Which kind is reported for a given invalid +input may change in any minor release, and new kinds may be added. +See [`PairwiseMappingFormat.ParserException`](@ref). + # Examples ```jldoctest julia> print(PAF.Errors.InvalidZero) @@ -392,8 +396,20 @@ using .Errors: Err The exception type thrown by iterating [`PAFReader`](@ref), or a failed `parse` of a [`PAFRecord`](@ref). The functions `try_parse` and `try_next!` may return (not throw) values of this type. -The line the error occurred may be accessed with the `.line` field. -The error kind may be accessed with the `.kind` field. + +The line the error occurred on may be accessed with the `.line` property. +For `try_parse` and `parse`, this is always 1. + +The `.kind` property is a [`PairwiseMappingFormat.Errors.Err`](@ref PairwiseMappingFormat.Errors) +describing why the line was rejected. It is intended for diagnostics, such as logging or +counting errors: An invalid line is often invalid in several ways, and which of them is +reported depends on implementation details. Therefore, the kind reported for a given +invalid input is deterministic within a version of this package, but may change in any +minor release, and new kinds may be added. +In contrast, whether input is rejected, and the line the error is reported on, are stable. + +The error message printed by `showerror` includes the approximate position of the error +within the line. This position is not part of the API. # Examples ```jldoctest @@ -763,7 +779,7 @@ end Base.IteratorSize(::Type{<:PAFReader}) = Base.SizeUnknown() Base.eltype(::Type{<:PAFReader}) = PAFRecord -Base.close(reader::PAFReader) = close(reader.io) # TODO: Docstring +Base.close(reader::PAFReader) = close(reader.io) function PAFReader(f::Function, io::Union{IO, AbstractBufReader}; kwargs...) reader = PAFReader(io; kwargs...) From 877b4c52f44ba6df53b014ceb39af02aefeb4f43 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 14:03:32 +0200 Subject: [PATCH 15/20] Allow biologically senseless, but spec-compliant records --- CHANGELOG.md | 4 ++++ docs/src/index.md | 7 +++++-- src/PairwiseMappingFormat.jl | 19 +++++++++---------- test/runtests.jl | 31 ++++++++++++++++++++++++++++++- 4 files changed, 48 insertions(+), 13 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 459f5eb..ea3f20c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -63,6 +63,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 skips all consecutive whitespace-only lines, and the error is reported at the first of them, instead of at the first non-whitespace byte after them. +- Empty query and target intervals, i.e. where the start equals the end, are now + accepted anywhere within the sequence. Previously, they were rejected with + `InvalidZero` at the start of the sequence, where the end is zero, but accepted elsewhere. + ### Removed - **Breaking:** The `PAFRecord(buffer_size::Int)` constructor. diff --git a/docs/src/index.md b/docs/src/index.md index 0a82629..ed4ec7f 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -105,11 +105,14 @@ underlying fields as they are stored in the `PAFRecord`. in the query and target, respectively. This does not take strandedness into account. * `qend` and `tend`. Of type `Int`. Rightmost (i.e. highest) position of the alignment in the query and target, respectively. This does not take strandedness into account. - As such, we always have `record.qend - record.qstart >= 0`, and likewise for the target. + The interval may be empty, in which case `record.qstart:record.qend` is an empty range. + As such, we always have `record.qend - record.qstart >= -1`, and likewise for the target. * `alnlen::Int`. Length of alignment. That includes gaps in query and target, and may therefore not match either `qend - qstart + 1`, or that of the target. * `matches::Int`. Number of residue matches (not mismatches) in the alignment. - `matches / alnlen` is the alignment identity and is between 0 and 1. + `matches / alnlen` is the alignment identity, which is between 0 and 1 for well-formed + records. The parser does not check that `matches <= alnlen`, or that `alnlen > 0`. + See [`aln_identity`](@ref). * `mapq::Union{Int, Nothing}`. Mapping quality, from 0:254, where 254 is the better quality. When calibrated, the probability the mapping is correct should be 10^(-mapq / 10). A value of `nothing` indicates the mapping quality is unavailable. diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 81d59c5..518540e 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -63,8 +63,9 @@ Use [`is_mapped`](@ref) to check if a record is mapped. ending positions of the alignments on the query and target strand. These uses one-based, closed interval indexing as usual in Julia, but unlike the underlying PAF format. - The ending positions are always >= the starting ones, and the ending positions - are always <= the query/target lengths. + The intervals may be empty, in which case the ending position is the starting + position minus one, so `qstart:qend` is an empty range. Hence, the ending positions + are always >= the starting ones minus one, and always <= the query/target lengths. * `matches::Int` gives the number of nucleotides / residues which match (i.e. are equal) in the alignment. * `alnlen::Int` gives the length of the alignment @@ -95,7 +96,7 @@ struct PAFRecord tstart::UInt32 tend::UInt32 # always >= tstart matches::UInt32 - alnlen::UInt32 # always >= 1 + alnlen::UInt32 # may be zero, and may be smaller than matches mapq::UInt8 # stored as 0xff for missing strand::UInt8 # 0x00: *, 0x01: +, 0x02: - # 16 bits of padding, which we can use for something else later @@ -272,7 +273,9 @@ end Return the nucleotide / residue identity, defined as the number of matches divided by the alignment length. -This is a number between 0 and 1. +For well-formed records, this is a number between 0 and 1. However, the parser does not +check that the number of matches is at most the alignment length, so this may be larger +than 1. If the alignment length is zero, this returns `NaN` or `Inf`. Throws an [`UnmappedError`](@ref) if the record is unmapped. # Examples @@ -563,13 +566,9 @@ function parse_line(mem::ImmutableMemoryView{UInt8})::Union{PAFRecord, ParserExc # The parse_x functions will return the next index to begin from (after the tab) (qname, i) = @? parse_str(mem, 1) (qlen, i) = @? parse_int(mem, i, false) + # Intervals may be empty, i.e. start == end, so the end may be zero (qstart, i) = @? parse_int(mem, i, true) - # qend may only be zero for unmapped records, e.g. minimap2's `--paf-no-hit` - # lines, so peek at the strand to check that before parsing it below. (qend, i) = @? parse_int(mem, i, true) - if iszero(qend) && (i > lastindex(mem) || @inbounds(mem[i]) != UInt8('*')) - return ParserException((i - 1) % Int32, Errors.InvalidZero) - end qend ≥ qstart || return ParserException(i % Int32, Errors.BackwardsIndices) qend > qlen && return ParserException(i % Int32, Errors.PositionOutOfBounds) @@ -595,7 +594,7 @@ function parse_line(mem::ImmutableMemoryView{UInt8})::Union{PAFRecord, ParserExc (tname, i) = @? parse_str(mem, i + 2) # i + 2 because strand plus tab took two bytes (tlen, i) = @? parse_int(mem, i, false) (tstart, i) = @? parse_int(mem, i, true) - (tend, i) = @? parse_int(mem, i, false) + (tend, i) = @? parse_int(mem, i, true) tend ≥ tstart || return ParserException(i % Int32, Errors.BackwardsIndices) tend > tlen && return ParserException(i % Int32, Errors.PositionOutOfBounds) (matches, i) = @? parse_int(mem, i, true) diff --git a/test/runtests.jl b/test/runtests.jl index c14b6ae..e934546 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -195,10 +195,14 @@ end @test with_replaced(good, allow_zero, "0") isa PAFRecord end - for not_allow_zero in [2, 4, 7, 9] + for not_allow_zero in [2, 7] @test with_replaced(good, not_allow_zero, "0") == Errors.InvalidZero end + # Ends may be zero, but only if the start is too, giving an empty interval + @test with_replaced(good, 4, "0") == Errors.BackwardsIndices + @test with_replaced(good, 9, "0") == Errors.BackwardsIndices + # qend and tend must be >= qstart tstart @test with_replaced(good, 3, "1000000") == Errors.BackwardsIndices @test with_replaced(good, 4, "2") == Errors.BackwardsIndices @@ -251,6 +255,31 @@ end @test (rec.qstart, rec.qend) == (Int(m) + 1, m) end +@testset "Empty intervals" begin + # Empty intervals are accepted anywhere within the sequence, including at its + # start and end. They become empty one-based ranges. + for (start, stop) in [(0, 0), (5, 5), (10, 10)] + rec = parse(PAFRecord, "q\t10\t$start\t$stop\t+\tt\t10\t$start\t$stop\t0\t0\t0") + @test (rec.qstart, rec.qend) == (start + 1, stop) + @test (rec.tstart, rec.tend) == (start + 1, stop) + @test isempty(rec.qstart:rec.qend) + @test isempty("0123456789"[rec.tstart:rec.tend]) + @test query_coverage(rec) == 0.0 + @test target_coverage(rec) == 0.0 + @test isnan(aln_identity(rec)) + end + + # But not out of bounds + @test try_parse("q\t10\t11\t11\t+\tt\t10\t0\t0\t0\t0\t0").kind == Errors.PositionOutOfBounds + @test try_parse("q\t10\t0\t0\t+\tt\t10\t11\t11\t0\t0\t0").kind == Errors.PositionOutOfBounds + + # matches and alnlen are not validated against each other + rec = parse(PAFRecord, "q\t10\t0\t5\t+\tt\t10\t0\t5\t7\t3\t0") + @test aln_identity(rec) == 7 / 3 + rec = parse(PAFRecord, "q\t10\t0\t5\t+\tt\t10\t0\t5\t2\t0\t0") + @test aln_identity(rec) == Inf +end + @testset "Bad reader" begin # Iterating a reader with bad data will throw good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" From 6196dcf4a0ac164827c327d24b58d09ce82f12d3 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 14:10:25 +0200 Subject: [PATCH 16/20] Add string, String and MemoryView constructors from record --- src/PairwiseMappingFormat.jl | 17 ++++++++++++++++- test/runtests.jl | 24 ++++++++++++++++++++++++ 2 files changed, 40 insertions(+), 1 deletion(-) diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 518540e..8c3390a 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -22,7 +22,7 @@ using PrecompileTools: @compile_workload, @setup_workload, @recompile_invalidati using StringViews: StringView end -using MemoryViews: ImmutableMemoryView +using MemoryViews: MemoryViews, ImmutableMemoryView using XAMAuxData.SAM: Auxiliary using BufferIO: AbstractBufReader, BufReader, ExpandReader, get_buffer, fill_buffer, consume @@ -46,6 +46,9 @@ julia> record.qname julia> record.qlen 301156 + +julia> startswith(String(record), "query1\t301156\t") +true ``` See also: [`PAFReader`](@ref) @@ -73,6 +76,13 @@ Use [`is_mapped`](@ref) to check if a record is mapped. this information is unavailable. * `is_rc::Bool` tells if the query and target strands are reverse-complement relative to each other. + +The line the record was parsed from, without the trailing newline, can be obtained +with `MemoryView(record)`, which returns an `ImmutableMemoryView{UInt8}` view of the +record's data, or `String(record)`, which copies it into a new `String`. +Unlike `MemoryView`, `String` and `string` do not alias the record's memory, so they +remain valid if the underlying memory is mutated. +For a non-copying string, use `StringView(MemoryView(record))`. """ struct PAFRecord # The line the record was parsed from, without the trailing newline. @@ -262,6 +272,11 @@ function tname(record::PAFRecord) return StringView(getfield(record, :data)[span]) end +# The original line, without the trailing newline +MemoryViews.MemoryView(record::PAFRecord) = getfield(record, :data) +Base.String(record::PAFRecord) = String(getfield(record, :data)) +Base.print(io::IO, record::PAFRecord) = (write(io, getfield(record, :data)); nothing) + # Note: We return nothing instead of missing because fuck missing function mapq(record::PAFRecord)::Union{Int, Nothing} x = getfield(record, :mapq) diff --git a/test/runtests.jl b/test/runtests.jl index e934546..dfd7b28 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -81,6 +81,30 @@ end end end +@testset "Raw line" begin + line = "q\t10\t0\t5\t+\tt\t20\t5\t10\t5\t5\t60\tab:Z:x" + for input in (line * "\n", line * "\r\n", line) + rec = parse(PAFRecord, input) + @test MemoryView(rec) isa ImmutableMemoryView{UInt8} + @test MemoryView(rec) == codeunits(line) + @test String(rec) isa String + @test String(rec) == line + @test string(rec) isa String + @test string(rec) == line + @test sprint(print, rec) == line + end + unmapped = "read1\t1523\t0\t0\t*\t*\t0\t0\t0\t0\t0\t0" + @test String(parse(PAFRecord, unmapped)) == unmapped + + # String and string copy the data, MemoryView does not + data = collect(codeunits(line)) + rec = parse(PAFRecord, data) + (s1, s2, mem) = (String(rec), string(rec), MemoryView(rec)) + data[1] = UInt8('Q') + @test s1 == s2 == line + @test first(mem) == UInt8('Q') +end + @testset "Coverage and identity" begin function get_rec(qlen::Int, tlen::Int, matches::Int, alnlen::Int) fields = split( From d148b19edfb56b302d503a057696ba673cdf19d5 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 14:17:47 +0200 Subject: [PATCH 17/20] Implement hash and == of records --- src/PairwiseMappingFormat.jl | 10 ++++++++++ test/runtests.jl | 24 ++++++++++++++++++++++++ 2 files changed, 34 insertions(+) diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 8c3390a..bd054c1 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -277,6 +277,16 @@ MemoryViews.MemoryView(record::PAFRecord) = getfield(record, :data) Base.String(record::PAFRecord) = String(getfield(record, :data)) Base.print(io::IO, record::PAFRecord) = (write(io, getfield(record, :data)); nothing) +# Records are equal if their lines are byte-equal, regardless of whether they share memory. +Base.:(==)(a::PAFRecord, b::PAFRecord) = getfield(a, :data) == getfield(b, :data) + +# Salt the hash, such that a record does not hash like its line as a string. +# Hash through StringView, because it is faster than the generic array hash. +const RECORD_HASH_SALT = 0x8a2e5d1f6b3c9047 % UInt +function Base.hash(record::PAFRecord, h::UInt) + return hash(StringView(getfield(record, :data)), h ⊻ RECORD_HASH_SALT) +end + # Note: We return nothing instead of missing because fuck missing function mapq(record::PAFRecord)::Union{Int, Nothing} x = getfield(record, :mapq) diff --git a/test/runtests.jl b/test/runtests.jl index dfd7b28..a625d02 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -105,6 +105,30 @@ end @test first(mem) == UInt8('Q') end +@testset "Equality and hashing" begin + line = "q\t10\t0\t5\t+\tt\t20\t5\t10\t5\t5\t60\tab:Z:x" + a = parse(PAFRecord, line) + for b in ( + parse(PAFRecord, line * "\n"), + parse(PAFRecord, line * "\r\n"), + parse(PAFRecord, collect(codeunits(line))), + copy(a), + ) + @test a == b + @test isequal(a, b) + @test hash(a) == hash(b) + end + c = parse(PAFRecord, line * "\tcd:i:1") + @test a != c + @test hash(a) != hash(c) + @test hash(a) != hash(line) + @test length(Set([a, copy(a), c])) == 2 + + unmapped = "read1\t1523\t0\t0\t*\t*\t0\t0\t0\t0\t0\t0" + @test parse(PAFRecord, unmapped) == parse(PAFRecord, unmapped * "\n") + @test parse(PAFRecord, unmapped) != a +end + @testset "Coverage and identity" begin function get_rec(qlen::Int, tlen::Int, matches::Int, alnlen::Int) fields = split( From 67c2f30d8727ab60084ca40ce37145d5f3f69b14 Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 14:24:11 +0200 Subject: [PATCH 18/20] Implement Base.peek(::PAFReader) --- docs/src/index.md | 1 + src/PairwiseMappingFormat.jl | 47 ++++++++++++++++++++++++++++++++++ test/runtests.jl | 49 ++++++++++++++++++++++++++++++++++++ 3 files changed, 97 insertions(+) diff --git a/docs/src/index.md b/docs/src/index.md index ed4ec7f..f0e8c6d 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -301,6 +301,7 @@ query_coverage target_coverage try_next! Base.isempty(::PAFReader) +Base.peek(::PAFReader) is_mapped PairwiseMappingFormat.Errors PairwiseMappingFormat.try_parse diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index bd054c1..696058c 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -886,7 +886,54 @@ true ``` """ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} + return expose(reader, next_view!(reader)) +end + +""" + peek(reader::PAFReader) -> Union{Nothing, PAFRecord, ParserException} + +Return what [`try_next!`](@ref) would return, without advancing the reader. +That is, return `nothing` if the reader is empty, a `ParserException` if the next record +is invalid, and the next record if it is valid. +Like `try_next!`, this does not throw on invalid records. + +While this does not advance the reader, if the reader was constructed with `copy=false`, +it invalidates any record previously returned by the reader, since it may need to +fill the reader's buffer. +A record returned by `peek` for a reader with `copy=false` is valid until the reader is +advanced, or `peek` or `isempty` is called again. +With `copy=true`, peeked records own their data, and are never invalidated. + +# Examples +```jldoctest +julia> reader = PAFReader(open(path_to_paf)); + +julia> peek(reader).qname +"query1" + +julia> try_next!(reader).qname +"query1" + +julia> peek(reader).qname +"another query" + +julia> close(reader) +``` +""" +function Base.peek(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} + # Invalidates the live view, if any, such that the buffer may be filled + consume_previous!(reader) + line = reader.line res = next_view!(reader) + # Undo the advancement. Since the previous line was consumed above, the + # buffer begins at the peeked line, and nothing needs to be consumed. + reader.to_consume = 0 + reader.line = line + return expose(reader, res) +end + +# Return a parsed view to the user, either as a copy or as a live view +function expose(reader::PAFReader, res::Union{Nothing, PAFRecord, ParserException}) res isa PAFRecord || return res # A copy owns its data, so the view never reaches the user reader.copy && return copy(res) diff --git a/test/runtests.jl b/test/runtests.jl index a625d02..28a173d 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -552,6 +552,55 @@ end @test only(collect(reader)).qname == "my_qname" end +@testset "peek" begin + good = open(readline, joinpath(@__DIR__, "example.paf")) + other = replace(good, "query1" => "other") + data = good * "\nbad\n\n \n" * other * "\r\n\n \n" + results(r) = r isa PAFRecord ? String(r) : r isa PairwiseMappingFormat.ParserException ? (r.kind, r.line) : r + expected = PAFReader(IOBuffer(data)) do reader + out = [] + while true + res = try_next!(reader) + push!(out, results(res)) + res === nothing && return out + end + end + @test length(expected) == 5 # good, bad, empty non-tailing line, other, nothing + + for copy in (true, false), chunk in (1, 3, 1000), n_peeks in (0, 1, 2) + chunks = [data[i:min(end, i + chunk - 1)] for i in 1:chunk:lastindex(data)] + reader = PAFReader(ChunkReader(chunks); copy) + for exp in expected + for _ in 1:n_peeks + @test results(peek(reader)) == exp + end + @test results(try_next!(reader)) == exp + end + @test peek(reader) === nothing + @test isempty(reader) + end + + # A peeked record is valid until the reader is advanced, also with copy=false + PAFReader(IOBuffer(good * "\n" * other); copy = false) do reader + rec = peek(reader) + @test !isempty(reader) + @test Base.isdone(reader) === false + @test String(rec) == good + @test String(try_next!(reader)) == good + @test String(peek(reader)) == other + @test String(try_next!(reader)) == other + @test peek(reader) === nothing + end + + # With copy=true, peeked records own their data + PAFReader(IOBuffer(good * "\n" * other)) do reader + rec = peek(reader) + try_next!(reader) + try_next!(reader) + @test String(rec) == good + end +end + @testset "isempty" begin good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" From 93ceba122d0e8d1fe8dd85236dbd0f62841ce34b Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 14:47:56 +0200 Subject: [PATCH 19/20] Perf: SWAR int parsing --- src/PairwiseMappingFormat.jl | 61 +++++++++++++++++++++++++++++++----- test/runtests.jl | 31 ++++++++++++++++++ 2 files changed, 84 insertions(+), 8 deletions(-) diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index 696058c..a7e75ff 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -692,18 +692,61 @@ function parse_str( end end -# TODO: 40% of runtime is spent in this function. Microoptimize it +# Return the number of leading bytes of `x` (in memory order) which are ASCII digits, +# from 0 to 8, where `x` is 8 bytes loaded as little-endian. +function n_leading_digits(x::UInt64)::Int + # A byte b is a digit iff the high nibbles of both b and b + 6 are 3. + # Adding 6 may carry into the next byte if b >= 0xfa, but that only corrupts bytes + # after a non-digit, which are not counted. + t = (x & 0xf0f0f0f0f0f0f0f0) | (((x + 0x0606060606060606) & 0xf0f0f0f0f0f0f0f0) >> 4) + return trailing_zeros(t ⊻ 0x3333333333333333) >>> 3 +end + +# Compute the value of the first `nd` digits of `x`, with `nd` in 0:8, and `x` +# as in `n_leading_digits`. +function swar_digits(x::UInt64, nd::Int)::Int + # Subtracting may borrow from bytes after the digits, which are shifted out. + # Shifting moves the digits to the most significant bytes, with zeros in front. + x = (x - 0x3030303030303030) << ((8 - nd) << 3) + x = ((x * 10) + (x >> 8)) & 0x00ff00ff00ff00ff + x = ((x * 100) + (x >> 16)) & 0x0000ffff0000ffff + x = ((x * 10000) + (x >> 32)) & 0x00000000ffffffff + return x % Int +end + +# Load 8 bytes beginning at index i as little-endian. Index i + 7 must be in bounds. +function load_u64(v::ImmutableMemoryView{UInt8}, i::Int)::UInt64 + return GC.@preserve v ltoh(unsafe_load(Ptr{UInt64}(pointer(v, i)))) +end + function parse_int( v::ImmutableMemoryView{UInt8}, from::Int, allow_zero::Bool, at_end::Bool = false, )::Union{Tuple{Int, Int}, ParserException} - # Accumulate in an Int, which can't overflow because we check against - # typemax(UInt32) after every digit. The result therefore fits in an UInt32. + # Fast path: Parse integers of up to 7 digits from 8 bytes at once, avoiding a + # data-dependent loop. Longer integers, i.e. >= 10 Mbp, continue in the loop below. + # A similar second fast path for 8-10 digits was measured to be slower. n = 0 - i = from - for outer i in from:lastindex(v) + start = from + if from + 7 <= lastindex(v) + x = load_u64(v, from) + nd = n_leading_digits(x) + n = swar_digits(x, nd) + if nd < 8 + # The byte after the digits is in bounds, since there are 8 bytes + i = from + nd + @inbounds(v[i]) == UInt8('\t') || return ParserException(i % Int32, Errors.BadInteger) + iszero(nd) && return ParserException(i % Int32, Errors.EmptyInteger) + (!allow_zero & iszero(n)) && return ParserException(i % Int32, Errors.InvalidZero) + return (n, i + 1) + end + start = from + 8 + end + # Slow path. Accumulate in an Int, which can't overflow because we check against + # typemax(UInt32) after every digit. The result therefore fits in an UInt32. + for i in start:lastindex(v) b = @inbounds v[i] # If we hit a tab, this field is over. if b == UInt8('\t') @@ -724,9 +767,9 @@ function parse_int( return if !at_end ParserException(lastindex(v) % Int32, Errors.TooFewFields) elseif from > lastindex(v) - ParserException(i % Int32, Errors.EmptyInteger) + ParserException(from % Int32, Errors.EmptyInteger) else - (n, i + 1) + (n, lastindex(v) + 1) end end @@ -946,7 +989,9 @@ function next_view!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserExceptio consume_previous!(reader) (buffer, last_index, line_len) = buffer_line(reader.io) line = @inbounds buffer[1:last_index] - res = parse_line(line) + # Inlining avoids returning the record through a union in memory, and copying it. + # This was measured to make the reader about 8% faster. + res = @inline parse_line(line) # A whitespace-only line, which can never be a valid record, is the end of the # stream if the rest of it is whitespace. Checking only on errors keeps the # check off the hot path. diff --git a/test/runtests.jl b/test/runtests.jl index 28a173d..fe63176 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -303,6 +303,37 @@ end @test (rec.qstart, rec.qend) == (Int(m) + 1, m) end +@testset "Integer parsing" begin + # Integers are parsed 8 bytes at a time, unless they have more than 7 digits, + # or are close to the end of the line. Test both cases for all digit counts. + good = "my_qname\t301156\t3\t301145\t+\tmy_tname\t6701780\t2764860\t3066002\t301142\t301142\t0" + m = Int(typemax(UInt32)) + alnlen(x) = x isa PAFRecord ? x.alnlen : x + for n_digits in 1:11, n_zeros in 0:9 + s = "0"^n_zeros * first("9182736450918273645", n_digits) + n = parse(Int, s) + expected = n > m ? Errors.IntegerOverflow : n + # alnlen is close to the end of the line, except when there is aux data + @test alnlen(with_replaced(good, 11, s)) == expected + @test alnlen(with_replaced(good * "\tcm:i:1", 11, s)) == expected + # Same with a bad byte after the integer + @test with_replaced(good, 11, s * "/") == (n > m ? Errors.IntegerOverflow : Errors.BadInteger) + end + + # A bad byte at every position of an integer, including bytes which may carry + # into the next byte when adding to them. + for n_digits in 1:12, pos in 1:n_digits, bad in (0x2f, 0x3a, 0x20, 0x00, 0xfa, 0xff, UInt8('a')) + bytes = collect(codeunits(first("123456789123", n_digits))) + bytes[pos] = bad + s = String(bytes) + for field in (10, 11, 12) + # 1234567891 fits in an UInt32, so only 11 digits overflow before the bad byte + expected = pos > 11 ? Errors.IntegerOverflow : Errors.BadInteger + @test with_replaced(good, field, s) == expected + end + end +end + @testset "Empty intervals" begin # Empty intervals are accepted anywhere within the sequence, including at its # start and end. They become empty one-based ranges. From c35b1f2a68af9bca08530909e3f40d07c5a1248b Mon Sep 17 00:00:00 2001 From: Jakob Nybo Andersen Date: Tue, 20 Oct 2026 21:20:29 +0200 Subject: [PATCH 20/20] Update docstrings --- docs/src/index.md | 13 +++---------- src/PairwiseMappingFormat.jl | 36 ++++++++++++++---------------------- 2 files changed, 17 insertions(+), 32 deletions(-) diff --git a/docs/src/index.md b/docs/src/index.md index f0e8c6d..ec79e9a 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -49,13 +49,10 @@ julia> PAFReader(open(path_to_paf)) do reader The [`PAFReader`](@ref) constructor takes an optional keyword `copy`, which defaults to `true`. If it is `false`, the [`PAFRecord`](@ref)s returned by iterating the reader, or by [`try_next!`](@ref), are views into the reader's buffer, instead of owning a copy of their data. -Such a record is only valid until the reader is advanced, and may be invalidated by -calling [`isempty`](@ref Base.isempty(::PAFReader)) on the reader. +Such a record is only valid until the reader is advanced, or otherwise mutated, which +may happen by other non-advancing calls such as `peek` or `isempty` on the reader. This reduces allocations and copying and gives a slight speed increase, but may result in bugs if records from previous iterations are stored. -Using an invalidated record is not detected, and does not throw an error. Instead, its data -may be silently overwritten with the content of later lines, so it may appear to be a -different, or even malformed, record. Use `copy` on a record to obtain a record that owns its data: ```jldoctest @@ -65,14 +62,10 @@ julia> PAFReader(open(path_to_paf); copy=false) do reader 2 ``` -!!! note - On Julia versions before 1.14, returning a `PAFRecord` from a function that is not - inlined allocates, even with `copy=false`. - At the moments, readers do not support seeking. ## Record -The immutable [`PAFRecord`](@ref) object represents a single line in a PAF file. +The [`PAFRecord`](@ref) object represents a single line in a PAF file. The individual columns of the PAF line is obtained by accessing properties of the records: ```jldoctest record diff --git a/src/PairwiseMappingFormat.jl b/src/PairwiseMappingFormat.jl index a7e75ff..cf1760c 100644 --- a/src/PairwiseMappingFormat.jl +++ b/src/PairwiseMappingFormat.jl @@ -65,10 +65,7 @@ Use [`is_mapped`](@ref) to check if a record is mapped. * `qstart`, `qend`, `tstart` and `tend::Int`, and give the starting and ending positions of the alignments on the query and target strand. These uses one-based, closed interval indexing as usual in Julia, but unlike - the underlying PAF format. - The intervals may be empty, in which case the ending position is the starting - position minus one, so `qstart:qend` is an empty range. Hence, the ending positions - are always >= the starting ones minus one, and always <= the query/target lengths. + the underlying PAF format. The intervals `qstart:qend` or `tstart:tend` may be empty. * `matches::Int` gives the number of nucleotides / residues which match (i.e. are equal) in the alignment. * `alnlen::Int` gives the length of the alignment @@ -142,7 +139,7 @@ end end """ - is_mapped(x::PAFRecord) -> Bool + is_mapped(x::PAFRecord)::Bool Compute whether the `PAFRecord` is mapped. Only the properties `qname` and `qlen`, and the auxiliary data of an unmapped record @@ -294,7 +291,7 @@ function mapq(record::PAFRecord)::Union{Int, Nothing} end """ - aln_identity(rec::Record) -> Float64 + aln_identity(rec::Record)::Float64 Return the nucleotide / residue identity, defined as the number of matches divided by the alignment length. @@ -317,7 +314,7 @@ function aln_identity(record::PAFRecord)::Float64 end """ - query_coverage(rec::Record) -> Float64 + query_coverage(rec::Record)::Float64 Compute the fraction of the query covered by the alignment. Throws an [`UnmappedError`](@ref) if the record is unmapped. @@ -336,7 +333,7 @@ function query_coverage(record::PAFRecord)::Float64 end """ - target_coverage(rec::Record) -> Float64 + target_coverage(rec::Record)::Float64 Compute the fraction of the target covered by the alignment. Throws an [`UnmappedError`](@ref) if the record is unmapped. @@ -355,7 +352,7 @@ function target_coverage(record::PAFRecord)::Float64 end """ - aux_data(rec::Record) -> XAMAuxData.SAM.Auxiliary + aux_data(rec::Record)::XAMAuxData.SAM.Auxiliary Return a lazily evaluated and validated `SAM.Auxiliary`, which is an `AbstractDict` of the auxiliary fields at the end of the record. @@ -513,7 +510,7 @@ macro var"?"(expr) end """ - try_parse(x) -> Union{PAFRecord, ParserException} + try_parse(x)::Union{PAFRecord, ParserException} Try parsing `x` into a `PAFRecord`. `x` may be any type that implements `MemoryView`, such as `String` or `Vector{UInt8}`. @@ -889,7 +886,7 @@ function Base.iterate(reader::PAFReader, (::Nothing) = nothing) end """ - try_next!(reader::PAFReader) -> Union{Nothing, PAFRecord, ParserException} + try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} Try to read a record from the [`PAFReader`](@ref) and advance the reader. This may yield one of three options: * If the reader is empty, return `nothing` and do not advance the reader @@ -933,19 +930,14 @@ function try_next!(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException end """ - peek(reader::PAFReader) -> Union{Nothing, PAFRecord, ParserException} + peek(reader::PAFReader)::Union{Nothing, PAFRecord, ParserException} Return what [`try_next!`](@ref) would return, without advancing the reader. -That is, return `nothing` if the reader is empty, a `ParserException` if the next record -is invalid, and the next record if it is valid. -Like `try_next!`, this does not throw on invalid records. While this does not advance the reader, if the reader was constructed with `copy=false`, -it invalidates any record previously returned by the reader, since it may need to -fill the reader's buffer. -A record returned by `peek` for a reader with `copy=false` is valid until the reader is -advanced, or `peek` or `isempty` is called again. -With `copy=true`, peeked records own their data, and are never invalidated. +it invalidates any record previously returned by the reader. Similarly, records +returned by the function with a reade with `copy=false` are invalidated like records +produced by iterating the reader. # Examples ```jldoctest @@ -1106,8 +1098,8 @@ Check whether `reader` has no more records, i.e. whether [`try_next!`](@ref) wou `copy=false`, it may invalidate the previous record by mutating the reader's inner state. In contrast, `Base.isdone(reader)`, never invalidates a record, -and returns `missing` if it can't determine whether the reader is done without invalidating -previous records. +and returns `missing` if it can't determine whether the reader is done +without invalidating previous records. """ function Base.isempty(reader::PAFReader) done = Base.isdone(reader)