From 118e6d567a328fa5af5ea52bb39a887a49d0e338 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Wed, 29 Jan 2025 16:15:50 +1000
Subject: [PATCH 01/16] Functions to set point format in the header

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/LASDatasets.jl  |  2 +-
 src/header.jl       | 60 +++++++++++++++++++++++++++++++--------------
 src/parse_points.jl |  4 +--
 src/points.jl       |  1 +
 test/header.jl      | 40 ++++++++++++++++++++++++++++--
 5 files changed, 84 insertions(+), 23 deletions(-)

diff --git a/src/LASDatasets.jl b/src/LASDatasets.jl
index 6c6babc..07aa6fb 100644
--- a/src/LASDatasets.jl
+++ b/src/LASDatasets.jl
@@ -39,7 +39,7 @@ export scan_direction, user_data
 # header interface
 export LasHeader, is_standard_gps, is_wkt, file_source_id, global_encoding, las_version, system_id, software_id, creation_day_of_year, creation_year
 export header_size, point_data_offset, point_record_length, record_format, point_format,  number_of_points, number_of_vlrs, number_of_evlrs, evlr_start, spatial_info, scale, num_return_channels
-export set_las_version!, set_spatial_info!, set_point_data_offset!, set_point_format!, set_point_record_length!, set_point_record_count!, set_num_vlr!, set_num_evlr!, set_gps_week_time_bit!
+export set_las_version!, set_point_format!, set_spatial_info!, set_point_data_offset!, set_point_format!, set_point_record_length!, set_point_record_count!, set_num_vlr!, set_num_evlr!, set_gps_week_time_bit!
 export set_gps_standard_time_bit!, is_internal_waveform, is_external_waveform, unset_wkt_bit!
 export set_waveform_external_bit!, set_waveform_internal_bit!, set_synthetic_return_numbers_bit!, unset_synthetic_return_numbers_bit!
 export set_wkt_bit!, get_number_of_points_by_return, set_number_of_points_by_return!, waveform_record_start
diff --git a/src/header.jl b/src/header.jl
index cfba31f..814dfc2 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -495,28 +495,47 @@ num_return_channels(h::LasHeader) = las_version(h) ≥ v"1.4" ? 15 : 5
 """
     $(TYPEDSIGNATURES)
 
-Get the LAS version in a header `h`
+Set the LAS specification version in a header `h` to version `v`
 """
 function set_las_version!(h::LasHeader, v::VersionNumber)
-    if any([
-        (v ≤ v"1.1") && (get_point_format_id(point_format(h)) ≥ 2),
-        (v ≤ v"1.2") && (get_point_format_id(point_format(h)) ≥ 4),
-        (v ≤ v"1.3") && (get_point_format_id(point_format(h)) ≥ 5),
-        (v ≤ v"1.4") && (get_point_format_id(point_format(h)) ≥ 11),
-        v ≥ v"1.5"
-    ])
-        error("Incompatible LAS version $(v) with point format $(point_format(h))")
-    end
+    _point_format_version_consistent(v, point_format(h))
     old_version = deepcopy(las_version(h))
     h.las_version = v
     h.header_size = get_header_size_from_version(v)
     h.data_offset += (h.header_size - get_header_size_from_version(old_version))
-    if (v ≥ v"1.4") && (get_point_format_id(point_format(h)) ≤ 5)
-        h.record_count = UInt64(h.legacy_record_count)
-        h.point_return_count = ntuple(i -> i ≤ 5 ? h.legacy_point_return_count[i] : 0, 15)
+    set_point_record_count!(h, number_of_points(h))
+    set_number_of_points_by_return!(h, get_number_of_points_by_return(h))
+    return nothing
+end
+
+function _point_format_version_consistent(v::VersionNumber, ::Type{TPoint}) where {TPoint <: LasPoint}
+    point_format_id = get_point_format_id(TPoint)
+    if any([
+        (v ≤ v"1.1") && (point_format_id ≥ 2),
+        (v ≤ v"1.2") && (point_format_id ≥ 4),
+        (v ≤ v"1.3") && (point_format_id ≥ 5),
+        (v ≤ v"1.4") && (point_format_id ≥ 11),
+        v ≥ v"1.5"
+    ])
+        error("Incompatible LAS version $(v) with point format $(point_format_id)")
     end
 end
 
+"""
+    $(TYPEDSIGNATURES)
+
+Set the point format in a header `h` to a new value, `TPoint`
+"""
+function set_point_format!(h::LasHeader, ::Type{TPoint}) where {TPoint <: LasPoint}
+    v = las_version(h)
+    _point_format_version_consistent(v, TPoint)
+    old_format_id = h.data_format_id
+    h.data_format_id = get_point_format_id(TPoint)
+    h.data_record_length += (byte_size(TPoint) - byte_size(LasPoint{Int(old_format_id)}))
+    set_point_record_count!(h, number_of_points(h))
+    set_number_of_points_by_return!(h, get_number_of_points_by_return(h))
+end
+
 """
     $(TYPEDSIGNATURES)
 
@@ -554,13 +573,16 @@ end
 Set the number of points in a LAS file with a header `header`
 """
 function set_point_record_count!(header::LasHeader, num_points::Integer)
-    if las_version(header) == v"1.4"
+    if (las_version(header) == v"1.4")
         @assert num_points ≤ typemax(UInt64) "Can't have more than $(typemax(UInt64)) points for LAS v1.4"
-        header.record_count = UInt64(num_points)
     else
         @assert num_points ≤ typemax(UInt32) "Can't have more than $(typemax(UInt32)) points for LAS v1.0-1.3"
+    end
+    header.record_count = UInt64(num_points)
+    if get_point_format_id(point_format(header)) ≤ 5
         header.legacy_record_count = UInt32(num_points)
     end
+    return nothing
 end
 
 """
@@ -725,10 +747,12 @@ function set_number_of_points_by_return!(header::LasHeader, points_per_return::N
     return_limit = las_version(header) ≥ v"1.4" ? typemax(UInt64) : typemax(UInt32)
     @assert all(points_per_return .≤ return_limit) "Maximum allowed points per return count is $return_limit"
     @assert N == num_return_channels(header) "Number of returns $N doesn't match what's in header $(num_return_channels(header))"
-    if las_version(header) ≥ v"1.4"
-        header.point_return_count = points_per_return
+    # note - internally, we store the point return count up to 15 places, even if the version spec only needs 5 - saves having to redefine field types for header
+    header.point_return_count = ntuple(i -> i ≤ N ? points_per_return[i] : 0, 15)
+    if get_point_format_id(point_format(header)) ≤ 5
+        header.legacy_point_return_count = ntuple(i -> i ≤ N ? points_per_return[i] : 0, 5)
     else
-        header.legacy_point_return_count = points_per_return
+        header.legacy_point_return_count = (0, 0, 0, 0, 0)
     end
 end
 
diff --git a/src/parse_points.jl b/src/parse_points.jl
index 6965a20..0f048c2 100644
--- a/src/parse_points.jl
+++ b/src/parse_points.jl
@@ -6,8 +6,8 @@ $(METHODLIST)
 function laspoint end
 laspoint(::Type{TPoint}, p::TPoint, xyz::SpatialInfo) where TPoint = p
 
-byte_size(point::TPoint) where {N, TPoint <: LasPoint{N}} = byte_size(TPoint)
-byte_size(::Type{TPoint}) where {N, TPoint <: LasPoint{N}} = sum(sizeof.(eltype.(fieldtypes(TPoint))))
+
+byte_size(::Type{TPoint}) where {TPoint <: LasPoint} = sum(sizeof.(eltype.(fieldtypes(get_point_format(TPoint)))))
 
 byte_size(vector::Type{SVector{N,T}}) where {N,T} = sizeof(T) * N
 
diff --git a/src/points.jl b/src/points.jl
index abcbeec..7dbb022 100644
--- a/src/points.jl
+++ b/src/points.jl
@@ -227,6 +227,7 @@ Get the concrete point format struct from an abstract `LasPoint` type
 
 $(METHODLIST)
 """
+get_point_format(::Type{TPoint}) where {TPoint <: LasPoint} = TPoint
 get_point_format(::Type{LasPoint{0}}) = LasPoint0
 get_point_format(::Type{LasPoint{1}}) = LasPoint1
 get_point_format(::Type{LasPoint{2}}) = LasPoint2
diff --git a/test/header.jl b/test/header.jl
index 5c970aa..b7b3997 100644
--- a/test/header.jl
+++ b/test/header.jl
@@ -6,7 +6,9 @@
         creation_year = UInt16(2024),
         data_format_id = 0x01,
         data_record_length = UInt16(28),
-        record_count = UInt64(100)
+        record_count = UInt64(100),
+        point_return_count = UInt64.((100, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)),
+        legacy_point_return_count = UInt32.((100, 0, 0, 0, 0))
     )
 
     # make sure we can get and set properties ok
@@ -27,6 +29,41 @@
     @test number_of_points(header) == 100
     set_point_record_count!(header, 150)
     @test number_of_points(header) == 150
+    @test get_number_of_points_by_return(header) == (100, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
+
+    # modifying record type and las version
+    set_las_version!(header, v"1.3")
+    @test las_version(header) == v"1.3"
+    @test header_size(header) == 235
+    @test point_data_offset(header) == 236
+    @test point_format(header) == LasPoint1
+    set_point_format!(header, LasPoint0)
+    @test point_format(header) == LasPoint0
+    # should now be 20 bytes (for las point 0) plus 2 undocumented extra bytes
+    @test point_record_length(header) == 22
+    # the point return counts should match the legacy counts here
+    @test get_number_of_points_by_return(header) == (100, 0, 0, 0, 0)
+
+    # we shouldn't be able to set an invalid point format
+    @test_throws ErrorException set_point_format!(header, LasPoint10)
+
+    # and we should be able to set the version back again
+    set_las_version!(header, v"1.4")
+    @test las_version(header) == v"1.4"
+    @test header_size(header) == 375
+    @test point_data_offset(header) == 376
+
+    
+    # but if we change to a newer LAS point format, we should get the same return counts, but set the legacy counts internally to 0
+    set_point_format!(header, LasPoint7)
+    @test get_number_of_points_by_return(header) == (100, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
+    # and the legacy counts should be 0
+    @test header.legacy_point_return_count == (0, 0, 0, 0, 0)
+
+    set_point_format!(header, LasPoint1)
+    @test point_format(header) == LasPoint1
+    @test point_record_length(header) == 30
+
     
     @test number_of_vlrs(header) == 0
     set_num_vlr!(header, 2)
@@ -36,7 +73,6 @@
     @test number_of_evlrs(header) == 1
 
     @test num_return_channels(header) == 15
-    @test get_number_of_points_by_return(header) == Tuple(zeros(Int, 15))
     points_per_return = ntuple(i -> 10, 15)
     set_number_of_points_by_return!(header, points_per_return)
     @test get_number_of_points_by_return(header) == points_per_return

From c1802a34e9fcc09168d2d91ce90b161f84af6a67 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 10:33:18 +1000
Subject: [PATCH 02/16] Combine LAS and user columns internally as one table

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/dataset.jl      | 119 ++++++++++++++++++++++++--------------------
 src/header.jl       |   1 +
 src/parse_points.jl |   6 +--
 src/write.jl        |   5 +-
 test/dataset.jl     |   7 +--
 5 files changed, 76 insertions(+), 62 deletions(-)

diff --git a/src/dataset.jl b/src/dataset.jl
index 1fe8f34..28e68d0 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -5,27 +5,24 @@ A wrapper around a LAS dataset. Contains point cloud data in tabular format as w
 
 $(TYPEDFIELDS)
 """
-mutable struct LASDataset
+struct LASDataset
     """The header from the LAS file the points were extracted from"""
-    const header::LasHeader
+    header::LasHeader
     
-    """LAS point cloud data stored in a Tabular format for convenience"""
-    const pointcloud::Table
-
-    """Additional user data assigned to each point that aren't standard LAS fields"""
-    _user_data::Union{Missing, FlexTable}
+    """Point cloud data stored in a Tabular format for convenience"""
+    pointcloud::FlexTable
     
     """Collection of Variable Length Records from the LAS file"""
-    const vlrs::Vector{LasVariableLengthRecord}
+    vlrs::Vector{LasVariableLengthRecord}
     
     """Collection of Extended Variable Length Records from the LAS file"""
-    const evlrs::Vector{LasVariableLengthRecord}
+    evlrs::Vector{LasVariableLengthRecord}
 
     """Extra user bytes packed between the Header block and the first VLR of the source LAS file"""
-    const user_defined_bytes::Vector{UInt8}
+    user_defined_bytes::Vector{UInt8}
 
     """Unit conversion factors applied to each axis when the dataset is ingested. This is reversed when you save the dataset to keep header/coordinate system information consistent"""
-    const unit_conversion::SVector{3, Float64}
+    unit_conversion::SVector{3, Float64}
 
     function LASDataset(header::LasHeader,
                         pointcloud::Table,
@@ -33,7 +30,7 @@ mutable struct LASDataset
                         evlrs::Vector{<:LasVariableLengthRecord},
                         user_defined_bytes::Vector{UInt8},
                         unit_conversion::SVector{3, Float64} = NO_CONVERSION)
-
+        pointcloud = FlexTable(pointcloud)
         # do a few checks to make sure everything is consistent between the header and other data
         point_format_from_table = get_point_format(pointcloud)
         point_format_from_header = point_format(header)
@@ -59,19 +56,16 @@ mutable struct LASDataset
 
         # if points don't have an ID column assigned to them, add it in
         if :id ∉ columnnames(pointcloud)
-            pointcloud = Table(pointcloud, id = collect(1:length(pointcloud)))
+            pointcloud.id = collect(1:length(pointcloud))
         end
 
         # need to split out our "standard" las columns from user-specific ones
         cols = collect(columnnames(pointcloud))
         these_are_las_cols = cols .∈ Ref(RECOGNISED_LAS_COLUMNS)
-        las_cols = cols[these_are_las_cols]
         other_cols = cols[.!these_are_las_cols]
-        las_pc = Table(NamedTuple{ (las_cols...,) }( (map(col -> getproperty(pointcloud, col), las_cols)...,) ))
         
-        make_consistent_header!(header, las_pc, vlrs, evlrs, user_defined_bytes)
+        make_consistent_header!(header, pointcloud, vlrs, evlrs, user_defined_bytes)
 
-        user_pc = isempty(other_cols) ? missing : FlexTable(NamedTuple{ (other_cols...,) }( (map(col -> getproperty(pointcloud, col), other_cols)...,) ))
         for col ∈ other_cols
             # account for potentially having an undocumented entry - in this case, don't add an ExtraBytes VLR
             if col != :undocumented_bytes
@@ -121,7 +115,7 @@ mutable struct LASDataset
                 end
             end
         end
-        return new(header, las_pc, user_pc, Vector{LasVariableLengthRecord}(vlrs), Vector{LasVariableLengthRecord}(evlrs), user_defined_bytes, unit_conversion)
+        return new(header, pointcloud, Vector{LasVariableLengthRecord}(vlrs), Vector{LasVariableLengthRecord}(evlrs), user_defined_bytes, unit_conversion)
     end
 end
 
@@ -158,11 +152,7 @@ end
 Extract point cloud data as a Table from a `LASDataset` `las`
 """
 function get_pointcloud(las::LASDataset)
-    if ismissing(las._user_data)
-        return las.pointcloud
-    else
-        return Table(las.pointcloud, las._user_data)
-    end
+    return las.pointcloud
 end
 
 """
@@ -215,9 +205,11 @@ function Base.show(io::IO, las::LASDataset)
     println(io, "LAS Dataset")
     println(io, "\tNum Points: $(length(get_pointcloud(las)))")
     println(io, "\tPoint Format: $(point_format(get_header(las)))")
-    println(io, "\tPoint Channels: $(columnnames(las.pointcloud))")
-    if !ismissing(las._user_data)
-        println(io, "\tUser Fields: $(columnnames(las._user_data))")
+    all_cols = columnnames(las.pointcloud)
+    is_las = all_cols .∈ Ref(RECOGNISED_LAS_COLUMNS)
+    println(io, "\tPoint Channels: $(all_cols[is_las])")
+    if any(.!is_las)
+        println(io, "\tUser Fields: $(all_cols[.!is_las])")
     end
     println(io, "\tVLR's: $(length(get_vlrs(las)))")
     println(io, "\tEVLR's: $(length(get_evlrs(las)))")
@@ -335,38 +327,51 @@ Add a column with name `column` and set of `values` to a `las` dataset
 function add_column!(las::LASDataset, column::Symbol, values::AbstractVector{T}) where T
     @assert length(values) == length(las.pointcloud) "Column size $(length(values)) inconsistent with number of points $(length(las.pointcloud))"
     check_user_type(T)
-    if ismissing(las._user_data)
-        las._user_data = FlexTable(NamedTuple{ (column,) }( (values,) ))
-    else
-        # make sure if we're replacing a column we correctly update the header size
-        if column ∈ columnnames(las._user_data)
-            las.header.data_record_length -= sizeof(eltype(getproperty(las._user_data, column)))
+    pointcloud = get_pointcloud(las)
+    
+    is_user = column ∉ RECOGNISED_LAS_COLUMNS
+    if is_user 
+        # need to update our header information and VLRs to track this user column
+        if column ∈ columnnames(pointcloud)
+            las.header.data_record_length -= sizeof(eltype(getproperty(pointcloud, column)))
         end
-        Base.setproperty!(las._user_data, column, values)
-    end
-    las.header.data_record_length += sizeof(T)
-    vlrs = get_vlrs(las)
-    extra_bytes_vlrs = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
-    @assert length(extra_bytes_vlrs) ≤ 1 "Found $(length(extra_bytes_vlrs)) Extra Bytes VLRs when we can only have a max of 1"
-    if isempty(extra_bytes_vlrs)
-        extra_bytes_vlr = LasVariableLengthRecord(LAS_SPEC_USER_ID, ID_EXTRABYTES, "Extra Bytes Records", ExtraBytesCollection())
-        # make sure we add it to the dataset to account for offsets in the header etc.
-        add_vlr!(las, extra_bytes_vlr)
-    else
-        extra_bytes_vlr = extra_bytes_vlrs[1]
-    end
-    if T <: SVector
-        # user field arrays have to be saved as sequential extra bytes records with names of the form "column [i]" (zero indexing encouraged)
-        split_col_name = split_column_name(column, length(T))
-        for i ∈ 1:length(T)
-            add_extra_bytes!(las, split_col_name[i], eltype(T), extra_bytes_vlr)
+        las.header.data_record_length += sizeof(T)
+        vlrs = get_vlrs(las)
+        extra_bytes_vlrs = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
+        @assert length(extra_bytes_vlrs) ≤ 1 "Found $(length(extra_bytes_vlrs)) Extra Bytes VLRs when we can only have a max of 1"
+        if isempty(extra_bytes_vlrs)
+            extra_bytes_vlr = LasVariableLengthRecord(LAS_SPEC_USER_ID, ID_EXTRABYTES, "Extra Bytes Records", ExtraBytesCollection())
+            # make sure we add it to the dataset to account for offsets in the header etc.
+            add_vlr!(las, extra_bytes_vlr)
+        else
+            extra_bytes_vlr = extra_bytes_vlrs[1]
         end
-    else
-        add_extra_bytes!(las, column, T, extra_bytes_vlr)
+        if T <: SVector
+            # user field arrays have to be saved as sequential extra bytes records with names of the form "column [i]" (zero indexing encouraged)
+            split_col_name = split_column_name(column, length(T))
+            for i ∈ 1:length(T)
+                add_extra_bytes!(las, split_col_name[i], eltype(T), extra_bytes_vlr)
+            end
+        else
+            add_extra_bytes!(las, column, T, extra_bytes_vlr)
+        end
+        # make sure we keep our offset to our first EVLR consistent now we've crammed more data in
+        update_evlr_offset!(las)
+    elseif column ∉ has_columns(point_format(get_header(las)))
+        # we're adding a new LAS column, which will necessitate changing the point format (and possibly the version)
+        las_cols = filter(c -> c ∈ RECOGNISED_LAS_COLUMNS, columnnames(pointcloud))
+        push!(las_cols, column)
+        new_format = get_point_format(las_cols)
+        @warn "Adding column $(column) to LAS dataset requires changing the point format to $(new_format), be warned!"
+        
+        set_point_format!(las, new_format)
     end
-    # make sure we keep our offset to our first EVLR consistent now we've crammed more data in
-    update_evlr_offset!(las)
-    nothing
+
+    # now actually write the values to the column
+    Base.setproperty!(pointcloud, column, values)
+    
+    
+    return nothing
 end
 
 """
@@ -424,4 +429,8 @@ function add_extra_bytes!(las::LASDataset, col_name::Symbol, ::Type{T}, extra_by
     add_extra_bytes_to_collection!(get_data(extra_bytes_vlr), col_name, T)
     header = get_header(las)
     header.data_offset += sizeof(ExtraBytes)
+end
+
+function set_point_format!(las::LASDataset, ::Type{TPoint}) where {TPoint <: LasPoint}
+    set_point_format!(get_header(las), TPoint)
 end
\ No newline at end of file
diff --git a/src/header.jl b/src/header.jl
index 814dfc2..b6d0deb 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -534,6 +534,7 @@ function set_point_format!(h::LasHeader, ::Type{TPoint}) where {TPoint <: LasPoi
     h.data_record_length += (byte_size(TPoint) - byte_size(LasPoint{Int(old_format_id)}))
     set_point_record_count!(h, number_of_points(h))
     set_number_of_points_by_return!(h, get_number_of_points_by_return(h))
+    return nothing
 end
 
 """
diff --git a/src/parse_points.jl b/src/parse_points.jl
index 0f048c2..cd7d658 100644
--- a/src/parse_points.jl
+++ b/src/parse_points.jl
@@ -16,9 +16,9 @@ byte_size(vector::Type{SVector{N,T}}) where {N,T} = sizeof(T) * N
 
 Get the minimum point format that is compatible with the contents of a point cloud in a `table`
 """
-function get_point_format(table::AbstractVector{<:NamedTuple})
-    
-    columns = columnnames(table)
+get_point_format(table::AbstractVector{<:NamedTuple}) = get_point_format(columnnames(table))
+
+function get_point_format(columns::Union{AbstractVector{Symbol}, NTuple{N, Symbol}}) where N
 
     # these are the known formats - the ones without wave packets
     columns_per_format = map(n -> has_columns(get_point_format(LasPoint{n})), collect(0:9))
diff --git a/src/write.jl b/src/write.jl
index e555712..92d94dc 100644
--- a/src/write.jl
+++ b/src/write.jl
@@ -96,7 +96,10 @@ function write_las(io::IO, las::LASDataset)
     this_point_format = point_format(header)
     xyz = spatial_info(header)
 
-    user_fields = ismissing(las._user_data) ? () : filter(c -> c != :undocumented_bytes, columnnames(las._user_data))
+    cols = collect(columnnames(pc))
+    these_are_las_cols = cols .∈ Ref(RECOGNISED_LAS_COLUMNS)
+    other_cols = cols[.!these_are_las_cols]
+    user_fields = filter(c -> c != :undocumented_bytes, other_cols)
 
     write(io, header)
 
diff --git a/test/dataset.jl b/test/dataset.jl
index 4cdc98c..6e69514 100644
--- a/test/dataset.jl
+++ b/test/dataset.jl
@@ -58,9 +58,7 @@
     this_header = get_header(las)
     @test point_record_length(this_header) == LASDatasets.byte_size(LasPoint1) + 10
     # our user fields should be populated in the dataset
-    @test las._user_data isa FlexTable
-    @test length(las._user_data) == num_points
-    @test sort(collect(columnnames(las._user_data))) == [:other_thing, :thing]
+    @test sort(collect(columnnames(las.pointcloud))) == [:classification, :gps_time, :id, :other_thing, :position, :thing]
     this_pc = get_pointcloud(las)
     @test this_pc.thing == spicy_pc.thing
     @test this_pc.other_thing == spicy_pc.other_thing
@@ -183,4 +181,7 @@
     remove_vlr!(las, superseded_comment)
     @test number_of_evlrs(get_header(las)) == 1
     @test get_evlrs(las) == [new_commment]
+
+    # las = LASDataset(pc)
+    
 end
\ No newline at end of file

From 29c47e1fb5d18e0006a74ea9bc9e7078788b402c Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 11:55:57 +1000
Subject: [PATCH 03/16] Update point format and version if adding incompatible
 las column

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/dataset.jl  |  8 ++++----
 src/header.jl   | 11 +++++++++--
 test/dataset.jl |  8 +++++++-
 test/header.jl  | 12 +++---------
 4 files changed, 23 insertions(+), 16 deletions(-)

diff --git a/src/dataset.jl b/src/dataset.jl
index 28e68d0..0eb75c9 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -206,7 +206,7 @@ function Base.show(io::IO, las::LASDataset)
     println(io, "\tNum Points: $(length(get_pointcloud(las)))")
     println(io, "\tPoint Format: $(point_format(get_header(las)))")
     all_cols = columnnames(las.pointcloud)
-    is_las = all_cols .∈ Ref(RECOGNISED_LAS_COLUMNS)
+    is_las = collect(all_cols .∈ Ref(RECOGNISED_LAS_COLUMNS))
     println(io, "\tPoint Channels: $(all_cols[is_las])")
     if any(.!is_las)
         println(io, "\tUser Fields: $(all_cols[.!is_las])")
@@ -326,11 +326,11 @@ Add a column with name `column` and set of `values` to a `las` dataset
 """
 function add_column!(las::LASDataset, column::Symbol, values::AbstractVector{T}) where T
     @assert length(values) == length(las.pointcloud) "Column size $(length(values)) inconsistent with number of points $(length(las.pointcloud))"
-    check_user_type(T)
     pointcloud = get_pointcloud(las)
     
     is_user = column ∉ RECOGNISED_LAS_COLUMNS
     if is_user 
+        check_user_type(T)
         # need to update our header information and VLRs to track this user column
         if column ∈ columnnames(pointcloud)
             las.header.data_record_length -= sizeof(eltype(getproperty(pointcloud, column)))
@@ -359,10 +359,10 @@ function add_column!(las::LASDataset, column::Symbol, values::AbstractVector{T})
         update_evlr_offset!(las)
     elseif column ∉ has_columns(point_format(get_header(las)))
         # we're adding a new LAS column, which will necessitate changing the point format (and possibly the version)
-        las_cols = filter(c -> c ∈ RECOGNISED_LAS_COLUMNS, columnnames(pointcloud))
+        las_cols = filter(c -> c ∈ RECOGNISED_LAS_COLUMNS, collect(columnnames(pointcloud)))
         push!(las_cols, column)
         new_format = get_point_format(las_cols)
-        @warn "Adding column $(column) to LAS dataset requires changing the point format to $(new_format), be warned!"
+        @warn "Changing point format to $(new_format) to allow the inclusion of LAS column $(column)"
         
         set_point_format!(las, new_format)
     end
diff --git a/src/header.jl b/src/header.jl
index b6d0deb..a20b55a 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -528,8 +528,15 @@ Set the point format in a header `h` to a new value, `TPoint`
 """
 function set_point_format!(h::LasHeader, ::Type{TPoint}) where {TPoint <: LasPoint}
     v = las_version(h)
-    _point_format_version_consistent(v, TPoint)
-    old_format_id = h.data_format_id
+    minimal_required_version = lasversion_for_point(TPoint)
+    old_format = point_format(h)
+    old_format_id = get_point_format_id(old_format)
+    # make sure that the LAS version in the header is consistent with the point format we want - upgrade if necessary, but let the user know
+    if v < minimal_required_version
+        @warn "Updating LAS version from $(v) to $(minimal_required_version) to accomodate changing point format from $(old_format) to $(TPoint)"
+        set_las_version!(h, minimal_required_version)
+    end
+    _point_format_version_consistent(las_version(h), TPoint)
     h.data_format_id = get_point_format_id(TPoint)
     h.data_record_length += (byte_size(TPoint) - byte_size(LasPoint{Int(old_format_id)}))
     set_point_record_count!(h, number_of_points(h))
diff --git a/test/dataset.jl b/test/dataset.jl
index 6e69514..447bb41 100644
--- a/test/dataset.jl
+++ b/test/dataset.jl
@@ -182,6 +182,12 @@
     @test number_of_evlrs(get_header(las)) == 1
     @test get_evlrs(las) == [new_commment]
 
-    # las = LASDataset(pc)
+    # test modifying point formats and versions
+    las = LASDataset(pc)
+
+    # if we add a LAS column that isn't covered by the current format, the point format (and possibly LAS version) should be updated in the header
+    add_column!(las, :overlap, falses(length(pc)))
+    @test point_format(get_header(las)) == LasPoint6
+    @test las_version(get_header(las)) == v"1.4"
     
 end
\ No newline at end of file
diff --git a/test/header.jl b/test/header.jl
index b7b3997..5bc3c19 100644
--- a/test/header.jl
+++ b/test/header.jl
@@ -44,20 +44,14 @@
     # the point return counts should match the legacy counts here
     @test get_number_of_points_by_return(header) == (100, 0, 0, 0, 0)
 
-    # we shouldn't be able to set an invalid point format
-    @test_throws ErrorException set_point_format!(header, LasPoint10)
-
-    # and we should be able to set the version back again
-    set_las_version!(header, v"1.4")
+    # if we request setting to a higher point format than the LAS version in the header supports, the version should be updated to the minimal compatible one
+    set_point_format!(header, LasPoint7)
     @test las_version(header) == v"1.4"
     @test header_size(header) == 375
     @test point_data_offset(header) == 376
-
     
-    # but if we change to a newer LAS point format, we should get the same return counts, but set the legacy counts internally to 0
-    set_point_format!(header, LasPoint7)
+    # now, we should get the same return counts, but set the legacy counts internally to 0
     @test get_number_of_points_by_return(header) == (100, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
-    # and the legacy counts should be 0
     @test header.legacy_point_return_count == (0, 0, 0, 0, 0)
 
     set_point_format!(header, LasPoint1)

From 54987f92b09f0fe273e2ea42fce93b04565a970c Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 12:04:04 +1000
Subject: [PATCH 04/16] Plumb header functions through to las datasets

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/dataset.jl | 36 ++++++++++++++++++++++++++++++++++++
 1 file changed, 36 insertions(+)

diff --git a/src/dataset.jl b/src/dataset.jl
index 0eb75c9..67c07d6 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -433,4 +433,40 @@ end
 
 function set_point_format!(las::LASDataset, ::Type{TPoint}) where {TPoint <: LasPoint}
     set_point_format!(get_header(las), TPoint)
+end
+
+# plumb through header util functions to act on a LASDataset for convenience
+for header_func ∈ (
+    :las_version,
+    :file_source_id,
+    :global_encoding, 
+    :system_id,
+    :software_id,
+    :creation_day_of_year,
+    :creation_year,
+    :header_size,
+    :point_data_offset, 
+    :point_record_length,
+    :point_format,
+    :number_of_points,
+    :number_of_vlrs,
+    :number_of_evlrs,
+    :evlr_start,
+    :spatrial_info,
+    :scale,
+    :num_return_channels,
+    :is_standard_gps,
+    :is_wkt,
+    :set_gps_standard_time_bit!,
+    :is_internal_waveform,
+    :set_waveform_internal_bit!,
+    :set_waveform_external_bit!,
+    :unset_waveform_bits!,
+    :set_synthetic_return_numbers_bit!,
+    :unset_synthetic_return_numbers_bit!,
+    :set_wkt_bit!,
+    :get_number_of_points_by_return,
+    :waveform_record_start
+    )
+    @eval $header_func(las::LASDataset) = $header_func(get_header(las))
 end
\ No newline at end of file

From ede3a351d09cccd4d13382362ab173cbc1845475 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 14:51:29 +1000
Subject: [PATCH 05/16] Add function to add points to dataset

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/LASDatasets.jl |  2 +-
 src/dataset.jl     | 16 +++++++++++++++-
 src/header.jl      | 18 +++++++++++++-----
 test/dataset.jl    | 18 +++++++++++++++++-
 4 files changed, 46 insertions(+), 8 deletions(-)

diff --git a/src/LASDatasets.jl b/src/LASDatasets.jl
index 07aa6fb..6ba6f12 100644
--- a/src/LASDatasets.jl
+++ b/src/LASDatasets.jl
@@ -55,7 +55,7 @@ export get_classes, get_description, set_description!
 export @register_vlr_type, read_vlr_data, extract_vlr_type
 
 export LASDataset, get_header, get_pointcloud, get_vlrs, get_evlrs, get_user_defined_bytes, get_unit_conversion
-export add_column!, merge_column!, add_vlr!, remove_vlr!, set_superseded!
+export add_column!, add_points!, merge_column!, add_vlr!, remove_vlr!, set_superseded!
 
 # I/O methods
 export load_las, load_pointcloud, save_las, load_header, load_vlrs
diff --git a/src/dataset.jl b/src/dataset.jl
index 67c07d6..bb3a197 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -435,6 +435,20 @@ function set_point_format!(las::LASDataset, ::Type{TPoint}) where {TPoint <: Las
     set_point_format!(get_header(las), TPoint)
 end
 
+"""
+    $(TYPEDSIGNATURES)
+
+Add a collection of `points` to a `LASDataset`, `las`. Updates header information to ensure dataset consistency
+"""
+function add_points!(las::LASDataset, points::AbstractVector{<:NamedTuple})
+    pc = get_pointcloud(las)
+    @assert all(columnnames(points) .∈ Ref(columnnames(pc))) "Point fields being appended $(columnnames(points)) must match those present in the dataset, $(columnnames(pc))"
+    append!(pc, points)
+    # make sure we update the header info too!
+    _consolidate_point_header_info!(get_header(las), pc)
+    return nothing
+end
+
 # plumb through header util functions to act on a LASDataset for convenience
 for header_func ∈ (
     :las_version,
@@ -452,7 +466,7 @@ for header_func ∈ (
     :number_of_vlrs,
     :number_of_evlrs,
     :evlr_start,
-    :spatrial_info,
+    :spatial_info,
     :scale,
     :num_return_channels,
     :is_standard_gps,
diff --git a/src/header.jl b/src/header.jl
index a20b55a..288a74a 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -836,12 +836,8 @@ function make_consistent_header!(header::LasHeader,
     point_data_offset = header_size + vlr_size + length(user_defined_bytes)
     
     set_point_data_offset!(header, point_data_offset)
-    set_spatial_info!(header, get_spatial_info(pointcloud; scale = scale(header)))
     
-    set_point_record_count!(header, length(pointcloud))
-    returns = (:returnnumber ∈ columnnames(pointcloud)) ? pointcloud.returnnumber : ones(Int, length(pointcloud))
-    points_per_return = ntuple(r -> count(returns .== r), num_return_channels(header))
-    set_number_of_points_by_return!(header, points_per_return)
+    _consolidate_point_header_info!(header, pointcloud)
 
     if !isempty(vlrs)
         set_num_vlr!(header, length(vlrs))
@@ -865,4 +861,16 @@ function make_consistent_header!(header::LasHeader,
         # don't want waveform bits set for point formats that don't have waveform data
         unset_waveform_bits!(header)
     end
+
+    return nothing
+end
+
+function _consolidate_point_header_info!(header::LasHeader, pointcloud::AbstractVector{<:NamedTuple})
+    set_spatial_info!(header, get_spatial_info(pointcloud; scale = scale(header)))
+    
+    set_point_record_count!(header, length(pointcloud))
+    returns = (:returnnumber ∈ columnnames(pointcloud)) ? pointcloud.returnnumber : ones(Int, length(pointcloud))
+    points_per_return = ntuple(r -> count(returns .== r), num_return_channels(header))
+    set_number_of_points_by_return!(header, points_per_return)
+    return nothing
 end
\ No newline at end of file
diff --git a/test/dataset.jl b/test/dataset.jl
index 447bb41..0e9a3a6 100644
--- a/test/dataset.jl
+++ b/test/dataset.jl
@@ -184,10 +184,26 @@
 
     # test modifying point formats and versions
     las = LASDataset(pc)
-
     # if we add a LAS column that isn't covered by the current format, the point format (and possibly LAS version) should be updated in the header
     add_column!(las, :overlap, falses(length(pc)))
     @test point_format(get_header(las)) == LasPoint6
     @test las_version(get_header(las)) == v"1.4"
+
+    # we should be able to add new points in too
+    new_points = Table(
+        position = rand(SVector{3, Float64}, 10),
+        classification = rand(UInt8, 10),
+        gps_time = rand(10),
+        id = length(pc):(length(pc) + 10),
+        overlap = falses(10)
+    )
+    add_points!(las, new_points)
+    pointcloud = get_pointcloud(las)
+    # check the point contents to make sure our new points are there
+    @test length(pointcloud) == 110
+    @test pointcloud[101:end] == new_points
+    # also make sure our header information is correctly set
+    @test number_of_points(las) == 110
+    @test spatial_info(las) == LASDatasets.get_spatial_info(pointcloud)
     
 end
\ No newline at end of file

From 04c03bc2adbc84fa101afc087b3f212f8fdbc150 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 14:53:26 +1000
Subject: [PATCH 06/16] AbstractVectors for spatial info funcs

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/spatial_info.jl | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/spatial_info.jl b/src/spatial_info.jl
index 9c0b165..52cfbc1 100644
--- a/src/spatial_info.jl
+++ b/src/spatial_info.jl
@@ -68,7 +68,7 @@ Base.write(io::IO, info::SpatialInfo) = write_struct(io, info)
 
 const DEFAULT_SPATIAL_INFO = SpatialInfo(AxisInfo(POINT_SCALE, POINT_SCALE, POINT_SCALE), AxisInfo(0.0, 0.0, 0.0), AxisInfo(Range(Inf, -Inf), Range(Inf, -Inf), Range(Inf, -Inf)))
 
-function bounding_box(points::Vector{SVector{3, T}}) where {T <: Real}
+function bounding_box(points::AbstractVector{SVector{3, T}}) where {T <: Real}
     x_min = typemax(T)
     x_max = typemin(T)
     y_min = typemax(T)
@@ -103,7 +103,7 @@ function bounding_box(points::Vector{SVector{3, T}}) where {T <: Real}
     return (; xmin = x_min, ymin = y_min, zmin = z_min, xmax = x_max, ymax = y_max, zmax = z_max)
 end
 
-function get_spatial_info(points::Vector{SVector{3, T}}; scale::T = T(POINT_SCALE)) where {T <: Real}
+function get_spatial_info(points::AbstractVector{SVector{3, T}}; scale::T = T(POINT_SCALE)) where {T <: Real}
     bb = bounding_box(points)
 
     x_offset = determine_offset(bb.xmin, bb.xmax, scale)

From de754e25c7ea551afd2b904e5c75ff6a361dd778 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 16:53:45 +1000
Subject: [PATCH 07/16] Handle special case of adding synthetics column

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/dataset.jl | 3 +++
 1 file changed, 3 insertions(+)

diff --git a/src/dataset.jl b/src/dataset.jl
index bb3a197..2573118 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -365,6 +365,9 @@ function add_column!(las::LASDataset, column::Symbol, values::AbstractVector{T})
         @warn "Changing point format to $(new_format) to allow the inclusion of LAS column $(column)"
         
         set_point_format!(las, new_format)
+
+        # special case - if we're adding synthetic points, need to set the synthetic flag in the header
+        set_synthetic_return_numbers_bit!(las)
     end
 
     # now actually write the values to the column

From 9e4e348733f77ede3c4c8efc497def46c7a7ce35 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 16:56:13 +1000
Subject: [PATCH 08/16] Also handle synthetic flag in header consolidation

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/header.jl | 5 +++++
 1 file changed, 5 insertions(+)

diff --git a/src/header.jl b/src/header.jl
index 288a74a..030ea28 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -872,5 +872,10 @@ function _consolidate_point_header_info!(header::LasHeader, pointcloud::Abstract
     returns = (:returnnumber ∈ columnnames(pointcloud)) ? pointcloud.returnnumber : ones(Int, length(pointcloud))
     points_per_return = ntuple(r -> count(returns .== r), num_return_channels(header))
     set_number_of_points_by_return!(header, points_per_return)
+    if :synthetic ∈ columnnames(pointcloud)
+        set_synthetic_return_numbers_bit!(header)
+    else
+        unset_synthetic_return_numbers_bit!(header)
+    end
     return nothing
 end
\ No newline at end of file

From 4938a82bf381148fe7885b25c155343475fb4aaa Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 17:08:25 +1000
Subject: [PATCH 09/16] Make extract_vlr_type return 1 vlr or nothing

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/dataset.jl | 18 ++++++++----------
 src/read.jl    |  3 +--
 src/vlrs.jl    | 11 ++++++++---
 src/write.jl   |  6 +++---
 4 files changed, 20 insertions(+), 18 deletions(-)

diff --git a/src/dataset.jl b/src/dataset.jl
index 2573118..da4c222 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -76,15 +76,12 @@ struct LASDataset
 
                 # grab information about the existing ExtraBytes VLRs - need to see if we need to update them or not
                 extra_bytes_vlr = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
-                @assert length(extra_bytes_vlr) ≤ 1 "Found multiple Extra Bytes VLRs in LAS file!"
-                if isempty(extra_bytes_vlr)
+                if isnothing(extra_bytes_vlr)
                     extra_bytes_vlr = LasVariableLengthRecord(LAS_SPEC_USER_ID, ID_EXTRABYTES, "Extra Bytes", ExtraBytesCollection())
                     # make sure we add the VLR to our collection and update any header info
                     push!(vlrs, extra_bytes_vlr)
                     header.n_vlr += 1
                     header.data_offset += sizeof(extra_bytes_vlr)
-                else
-                    extra_bytes_vlr = extra_bytes_vlr[1]
                 end
                 extra_bytes_data = get_extra_bytes(get_data(extra_bytes_vlr))
                 user_field_names = Symbol.(name.(extra_bytes_data))
@@ -256,12 +253,16 @@ Add a `vlr` into the set of VLRs in a LAS dataset `las`.
 Note that this will modify the header content of `las`, including updating its LAS version to v1.4 if `vlr` is extended
 """
 function add_vlr!(las::LASDataset, vlr::LasVariableLengthRecord)
-    if is_extended(vlr) && isempty(get_evlrs(las))
+    if is_extended(vlr) && las_version(las) < v"1.4"
         # evlrs only supported in LAS 1.4
+        @warn "Upgrading LAS spec version to 1.4.0 from $(las_version(las)) to support use of EVLRs"
         set_las_version!(get_header(las), v"1.4")
     end
 
     header = get_header(las)
+
+    existing_vlr = extract_vlr_type(get_vlrs(las), get_user_id(vlr), get_record_id(vlr))
+    @assert isnothing(existing_vlr) "We already have a VLR with user ID $(get_user_id(vlr)) and record ID $(get_record_id(vlr)) in the las dataset!"
     
     if is_extended(vlr)
         set_num_evlr!(header, number_of_evlrs(header) + 1)
@@ -337,14 +338,11 @@ function add_column!(las::LASDataset, column::Symbol, values::AbstractVector{T})
         end
         las.header.data_record_length += sizeof(T)
         vlrs = get_vlrs(las)
-        extra_bytes_vlrs = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
-        @assert length(extra_bytes_vlrs) ≤ 1 "Found $(length(extra_bytes_vlrs)) Extra Bytes VLRs when we can only have a max of 1"
-        if isempty(extra_bytes_vlrs)
+        extra_bytes_vlr = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
+        if isnothing(extra_bytes_vlr)
             extra_bytes_vlr = LasVariableLengthRecord(LAS_SPEC_USER_ID, ID_EXTRABYTES, "Extra Bytes Records", ExtraBytesCollection())
             # make sure we add it to the dataset to account for offsets in the header etc.
             add_vlr!(las, extra_bytes_vlr)
-        else
-            extra_bytes_vlr = extra_bytes_vlrs[1]
         end
         if T <: SVector
             # user field arrays have to be saved as sequential extra bytes records with names of the form "column [i]" (zero indexing encouraged)
diff --git a/src/read.jl b/src/read.jl
index 6a89d25..fa1d5f4 100644
--- a/src/read.jl
+++ b/src/read.jl
@@ -107,8 +107,7 @@ function read_las_data(io::TIO, required_columns::TTuple=DEFAULT_LAS_COLUMNS;
     user_defined_bytes = read(io, header.data_offset - pos)
 
     extra_bytes_vlr = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
-    @assert length(extra_bytes_vlr) ≤ 1 "Found multiple extra bytes columns!"
-    extra_bytes = isempty(extra_bytes_vlr) ? ExtraBytes[] : get_extra_bytes(get_data(extra_bytes_vlr[1]))
+    extra_bytes = isnothing(extra_bytes_vlr) ? ExtraBytes[] : get_extra_bytes(get_data(extra_bytes_vlr))
 
     this_format = record_format(header, extra_bytes)
     xyz = spatial_info(header)
diff --git a/src/vlrs.jl b/src/vlrs.jl
index 770de04..2f9a4ad 100644
--- a/src/vlrs.jl
+++ b/src/vlrs.jl
@@ -259,9 +259,14 @@ is_srs(vlr::LasVariableLengthRecord) = vlr.record_id in (
 """
     $(TYPEDSIGNATURES)
 
-Extract all VLRs with a `user_id` and `record_id` from a collection of VLRs, `vlrs`
+Extract the VLR with a `user_id` and `record_id` from a collection of VLRs, `vlrs`
 """
 function extract_vlr_type(vlrs::Vector{<:LasVariableLengthRecord}, user_id::String, record_id::Integer)
-    matches_ids = findall(vlr -> (get_user_id(vlr) == user_id) && (get_record_id(vlr) == record_id), vlrs)
-    return vlrs[matches_ids]
+    # there should only ever be 1 VLR that has a particular user/record ID combination
+    matches_ids = findfirst(vlr -> (get_user_id(vlr) == user_id) && (get_record_id(vlr) == record_id), vlrs)
+    if isnothing(matches_ids)
+        return nothing
+    else
+        return vlrs[matches_ids]
+    end
 end
\ No newline at end of file
diff --git a/src/write.jl b/src/write.jl
index 92d94dc..a8f47c7 100644
--- a/src/write.jl
+++ b/src/write.jl
@@ -154,10 +154,10 @@ function get_record_bytes(records::StructVector{TRecord}, vlrs::Vector{LasVariab
 
     if user_field_bytes > 0
         # need to write the extra bytes fields in the same order as they appear in the VLR
-        extra_bytes_vlrs = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
-        @assert length(extra_bytes_vlrs) == 1 "Expected to find 1 Extra Bytes VLR, instead found $(length(extra_bytes_vlrs))"
+        extra_bytes_vlr = extract_vlr_type(vlrs, LAS_SPEC_USER_ID, ID_EXTRABYTES)
+        @assert !isnothing(extra_bytes_vlr) "Expected to find an Extra Bytes VLR but found none!"
         # get the order they appear in the VLR
-        user_field_names = unique(get_base_field_name.(Symbol.(name.(get_extra_bytes(get_data(extra_bytes_vlrs[1]))))))
+        user_field_names = unique(get_base_field_name.(Symbol.(name.(get_extra_bytes(get_data(extra_bytes_vlr))))))
         # create a mapping between the order in the VLR and the order in the record
         per_record_user_field_names = get_user_field_names(TRecord)
         user_field_idxs = indexin(user_field_names, collect(per_record_user_field_names))

From 3835bbe72a37b7a482ab2f355ba1f490608916da Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 17:13:05 +1000
Subject: [PATCH 10/16] Fix superseded case in tests

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 test/dataset.jl | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/test/dataset.jl b/test/dataset.jl
index 0e9a3a6..0ef4003 100644
--- a/test/dataset.jl
+++ b/test/dataset.jl
@@ -117,9 +117,11 @@
         "Text Area Description", 
         TextAreaDescription("This is the new dataset description")
     )
-    add_vlr!(las, new_desc)
     # mark the old one as superseded
     set_superseded!(las, desc)
+    # and add the new one
+    add_vlr!(las, new_desc)
+    
     vlrs = get_vlrs(las)
     @test length(vlrs) == 3
     @test vlrs[3] == new_desc

From ced6327834a41755b17ca19e670b3bd0980be749 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Thu, 30 Jan 2025 17:28:42 +1000
Subject: [PATCH 11/16] Correct functionality for scaling over all axes

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/header.jl       | 20 ++++++++++++--------
 src/spatial_info.jl | 30 ++++++++++++++++--------------
 2 files changed, 28 insertions(+), 22 deletions(-)

diff --git a/src/header.jl b/src/header.jl
index 030ea28..f847cdb 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -476,15 +476,19 @@ spatial_info(h::LasHeader) = h.spatial_info
 """
     $(TYPEDSIGNATURES)
 
-Get the scale for point positions in a LAS file from a header `h`. Checks consistent scale for ALL axes.
+Get the scale for point positions in a LAS file from a header `h` along an `axis` (x, y or z)
 """
-function scale(h::LasHeader) 
-
-    @assert (h.spatial_info.scale.x == h.spatial_info.scale.y) && (h.spatial_info.scale.y == h.spatial_info.scale.z) "We expect all axes to be scaled similarly"
-    
-    return h.spatial_info.scale.x
+function scale(h::LasHeader, axis::Symbol) 
+    return getproperty(h.spatial_info.scale, axis)
 end
 
+"""
+    $(TYPEDSIGNATURES)
+
+Get the scale for point positions in a LAS file from a header `h` along all axes
+"""
+scale(h::LasHeader) = h.spatial_info.scale
+
 """
     $(TYPEDSIGNATURES)
 
@@ -804,7 +808,7 @@ function make_consistent_header(pointcloud::AbstractVector{<:NamedTuple},
                                 vlrs::Vector{<:LasVariableLengthRecord}, 
                                 evlrs::Vector{<:LasVariableLengthRecord},
                                 user_defined_bytes::Vector{UInt8},
-                                scale::Real) where {TPoint <: LasPoint}
+                                scale::Union{Real, SVector{3, <:Real}}) where {TPoint <: LasPoint}
     version = lasversion_for_point(point_format)
 
     spatial_info = get_spatial_info(pointcloud; scale = scale)
@@ -866,7 +870,7 @@ function make_consistent_header!(header::LasHeader,
 end
 
 function _consolidate_point_header_info!(header::LasHeader, pointcloud::AbstractVector{<:NamedTuple})
-    set_spatial_info!(header, get_spatial_info(pointcloud; scale = scale(header)))
+    set_spatial_info!(header, get_spatial_info(pointcloud, scale(header)))
     
     set_point_record_count!(header, length(pointcloud))
     returns = (:returnnumber ∈ columnnames(pointcloud)) ? pointcloud.returnnumber : ones(Int, length(pointcloud))
diff --git a/src/spatial_info.jl b/src/spatial_info.jl
index 52cfbc1..adb57d3 100644
--- a/src/spatial_info.jl
+++ b/src/spatial_info.jl
@@ -103,28 +103,30 @@ function bounding_box(points::AbstractVector{SVector{3, T}}) where {T <: Real}
     return (; xmin = x_min, ymin = y_min, zmin = z_min, xmax = x_max, ymax = y_max, zmax = z_max)
 end
 
-function get_spatial_info(points::AbstractVector{SVector{3, T}}; scale::T = T(POINT_SCALE)) where {T <: Real}
+get_spatial_info(points::AbstractVector{SVector{3, T}}; scale::T = T(POINT_SCALE)) where {T <: Real} = get_spatial_info(points, AxisInfo{T}(scale, scale, scale))
+get_spatial_info(points::AbstractVector{SVector{3, T}}, scale::SVector{3, T}) where {T <: Real} = get_spatial_info(points, AxisInfo{T}(scale.x, scale.y, scale.z))
+
+function get_spatial_info(points::AbstractVector{SVector{3, T}}, scale::AxisInfo{T}) where {T <: Real}
     bb = bounding_box(points)
 
-    x_offset = determine_offset(bb.xmin, bb.xmax, scale)
-    y_offset = determine_offset(bb.ymin, bb.ymax, scale)
-    z_offset = determine_offset(bb.zmin, bb.zmax, scale)
+    x_offset = determine_offset(bb.xmin, bb.xmax, scale.x)
+    y_offset = determine_offset(bb.ymin, bb.ymax, scale.y)
+    z_offset = determine_offset(bb.zmin, bb.zmax, scale.z)
 
     offset = AxisInfo(x_offset, y_offset, z_offset)
 
-    x_min = scale * round((bb.xmin / scale) - 0.5)
-    y_min = scale * round((bb.ymin / scale) - 0.5)
-    z_min = scale * round((bb.zmin / scale) - 0.5)
-    x_max = scale * round((bb.xmax / scale) + 0.5)
-    y_max = scale * round((bb.ymax / scale) + 0.5)
-    z_max = scale * round((bb.zmax / scale) + 0.5)
-    
-    scale_info = AxisInfo(scale, scale, scale)   
+    x_min = scale.x * round((bb.xmin / scale.x) - 0.5)
+    y_min = scale.y * round((bb.ymin / scale.y) - 0.5)
+    z_min = scale.z * round((bb.zmin / scale.z) - 0.5)
+    x_max = scale.x * round((bb.xmax / scale.x) + 0.5)
+    y_max = scale.y * round((bb.ymax / scale.y) + 0.5)
+    z_max = scale.z * round((bb.zmax / scale.z) + 0.5) 
 
-    return SpatialInfo(scale_info, offset, AxisInfo(Range(x_max, x_min), Range(y_max, y_min), Range(z_max, z_min)))
+    return SpatialInfo(scale, offset, AxisInfo(Range(x_max, x_min), Range(y_max, y_min), Range(z_max, z_min)))
 end
 
-get_spatial_info(pc::AbstractVector{<:NamedTuple}; kwargs...) = get_spatial_info(pc.position; kwargs...)
+get_spatial_info(pc::Union{AbstractVector{<:NamedTuple}, FlexTable}; kwargs...) = get_spatial_info(pc.position; kwargs...)
+get_spatial_info(pc::Union{AbstractVector{<:NamedTuple}, FlexTable}, scale::Union{SVector, AxisInfo}) = get_spatial_info(pc.position, scale)
 
 function determine_offset(min_value, max_value, scale; threshold=10^7)
     s = round(Int64, ((min_value + max_value) / 2) / scale / threshold)

From ac06b7a0823604dc251a1ebf9b7384ddf63d574c Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Fri, 31 Jan 2025 09:16:53 +1000
Subject: [PATCH 12/16] Add docs + fix logic in add_points!

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 docs/src/header.md    |  1 +
 docs/src/interface.md | 14 ++++++++++++--
 docs/src/vlrs.md      |  4 ++--
 src/dataset.jl        |  8 ++++++--
 4 files changed, 21 insertions(+), 6 deletions(-)

diff --git a/docs/src/header.md b/docs/src/header.md
index 3bdace3..b841a41 100644
--- a/docs/src/header.md
+++ b/docs/src/header.md
@@ -38,6 +38,7 @@ You can also modify certain fields in the header, but one should note that for s
 
 ```@docs; canonical = false
 set_las_version!
+set_point_format!
 set_spatial_info!
 set_point_data_offset!
 set_point_record_length!
diff --git a/docs/src/interface.md b/docs/src/interface.md
index 95bd78f..cb0ad7a 100644
--- a/docs/src/interface.md
+++ b/docs/src/interface.md
@@ -70,18 +70,28 @@ Alternatively, if you just have the point cloud data as a Table:
 using StaticArrays
 using TypedTables
 
-pc = Table(position = rand(SVector{3, Float64}, 10), classification = rand(UIn8, 10))
+pc = Table(position = rand(SVector{3, Float64}, 10), classification = rand(UInt8, 10))
 save_las("my_las.las", pc)
 ```
 
 Note that when you supply just the point cloud outside of a `LASDataset`, *LASDatasets.jl* will automatically construct the appropriate header for you so you don't need to worry about the specifics of appropriate point formats etc. 
 
 ## Modifying LAS Contents
-You can modify point fields in your `LASDataset` by adding new columns or merging in values from an existing vector.
+You can modify point fields in your `LASDataset` by adding new columns or merging in values from an existing vector. Additionally, you can add new points to your dataset, however one should note that the user is responsible for correctly setting the appropriate fields such as synthetic flags themselves.
 
 ```@docs; canonical = false
 add_column!
 merge_column!
+add_points!
+```
+
+For example, if you want to add a set of synthetic points to your dataset, you can run:
+```julia
+las = load_las("my_las.las")
+# note - we need to set a synthetic column here for the existing points before we append points with this field
+add_column!(las, :synthetic, falses(number_of_points(las)))
+synthetic_points = Table(position = rand(SVector{3, Float64}, 5), classification = rand(UInt8, 5), synthetic = trues(5))
+add_points!(las, synthetic_points)
 ```
 
 You can also add or remove *(E)VLRs* using the following functions, and set an existing *(E)VLR* as *superseded* if it's an old copy of a record.
diff --git a/docs/src/vlrs.md b/docs/src/vlrs.md
index 8ece1bf..e050897 100644
--- a/docs/src/vlrs.md
+++ b/docs/src/vlrs.md
@@ -62,7 +62,7 @@ You can then read the *LAS* data and extract the classification lookup:
 las = load_las("pc.las")
 # look for the classification lookup VLR by checking for its user and record IDs
 lookup_vlr = extract_vlr_type(get_vlrs(las), "LASF_Spec", 0)
-lookup = get_data(lookup_vlr[1])
+lookup = get_data(lookup_vlr)
 ```
 
 ### Text Area Descriptions
@@ -150,7 +150,7 @@ So in our example, we can tell the system that records containing data of type `
 @register_vlr_type MyType "My Custom Records" collect(1:100)
 ```
 
-And now we can save our `MyType` *VLRs* into a *LAS* file in the same way as we did above for the register *VLR* types. Note that you can use the function `extract_vlr_type` on your collection of *VLRs* to pull out any *VLRs* with a specific user ID and record ID.
+And now we can save our `MyType` *VLRs* into a *LAS* file in the same way as we did above for the register *VLR* types. Note that you can use the function `extract_vlr_type` on your collection of *VLRs* to pull out the *VLR* with a specific user ID and record ID. 
 
 ```@docs; canonical = false
 extract_vlr_type
diff --git a/src/dataset.jl b/src/dataset.jl
index da4c222..2287c06 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -443,8 +443,12 @@ Add a collection of `points` to a `LASDataset`, `las`. Updates header informatio
 """
 function add_points!(las::LASDataset, points::AbstractVector{<:NamedTuple})
     pc = get_pointcloud(las)
-    @assert all(columnnames(points) .∈ Ref(columnnames(pc))) "Point fields being appended $(columnnames(points)) must match those present in the dataset, $(columnnames(pc))"
-    append!(pc, points)
+    # make new points a FlexTable so we can add any missing columns with 0 values to ensure consistency
+    missing_cols = filter(c -> c ∉ columnnames(points), columnnames(pc))
+    @warn "Adding default entries for missing columns $(missing_cols)"
+    # need to make sure columns are in the same order as in pc to avoid errors from TypedTables
+    new_points = Table(NamedTuple{ (columnnames(pc)...,) }( (map(col -> hasproperty(points, col) ? getproperty(points, col) : zeros(eltype(getproperty(pc, col)), length(points)), columnnames(pc))...,) ))
+    append!(pc, new_points)
     # make sure we update the header info too!
     _consolidate_point_header_info!(get_header(las), pc)
     return nothing

From 048126a31fcb4412f630a4a42d5c6a16f27e245b Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Fri, 31 Jan 2025 09:17:52 +1000
Subject: [PATCH 13/16] PR Feedback

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/dataset.jl | 1 -
 src/header.jl  | 1 +
 2 files changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/dataset.jl b/src/dataset.jl
index 2287c06..ec92594 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -371,7 +371,6 @@ function add_column!(las::LASDataset, column::Symbol, values::AbstractVector{T})
     # now actually write the values to the column
     Base.setproperty!(pointcloud, column, values)
     
-    
     return nothing
 end
 
diff --git a/src/header.jl b/src/header.jl
index f847cdb..694e616 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -766,6 +766,7 @@ function set_number_of_points_by_return!(header::LasHeader, points_per_return::N
     else
         header.legacy_point_return_count = (0, 0, 0, 0, 0)
     end
+    return nothing
 end
 
 """

From 92cf95db46f28bb77a8a5609c762d156e2a4e308 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Fri, 31 Jan 2025 09:24:00 +1000
Subject: [PATCH 14/16] Version bump to 0.4.0

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 Project.toml    | 2 +-
 test/dataset.jl | 7 +++++--
 2 files changed, 6 insertions(+), 3 deletions(-)

diff --git a/Project.toml b/Project.toml
index 307ea24..bf77b3c 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "LASDatasets"
 uuid = "cc498e2a-d443-4943-8f26-2a8a0f3c7cdb"
 authors = ["BenCurran98 <b.curran@fugro.com>"]
-version = "0.3.3"
+version = "0.4.0"
 
 [deps]
 BufferedStreams = "e1450e63-4bb3-523b-b2a4-4ffa8c0fd77d"
diff --git a/test/dataset.jl b/test/dataset.jl
index 0ef4003..46e8abd 100644
--- a/test/dataset.jl
+++ b/test/dataset.jl
@@ -193,17 +193,20 @@
 
     # we should be able to add new points in too
     new_points = Table(
+        id = length(pc):(length(pc) + 10),
         position = rand(SVector{3, Float64}, 10),
         classification = rand(UInt8, 10),
         gps_time = rand(10),
-        id = length(pc):(length(pc) + 10),
         overlap = falses(10)
     )
     add_points!(las, new_points)
     pointcloud = get_pointcloud(las)
     # check the point contents to make sure our new points are there
     @test length(pointcloud) == 110
-    @test pointcloud[101:end] == new_points
+    # equality checks are annoying when the column order gets switched internally, so check each column individually
+    for col ∈ columnnames(new_points)
+        @test getproperty(pointcloud, col)[101:end] == getproperty(new_points, col)
+    end
     # also make sure our header information is correctly set
     @test number_of_points(las) == 110
     @test spatial_info(las) == LASDatasets.get_spatial_info(pointcloud)

From 560dee9fac65ab7cab0e2b57faef5019b1a791f0 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Fri, 31 Jan 2025 10:15:55 +1000
Subject: [PATCH 15/16] Only log warning if there are missing cols

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 src/dataset.jl | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/src/dataset.jl b/src/dataset.jl
index ec92594..abb59cf 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -444,7 +444,9 @@ function add_points!(las::LASDataset, points::AbstractVector{<:NamedTuple})
     pc = get_pointcloud(las)
     # make new points a FlexTable so we can add any missing columns with 0 values to ensure consistency
     missing_cols = filter(c -> c ∉ columnnames(points), columnnames(pc))
-    @warn "Adding default entries for missing columns $(missing_cols)"
+    if !isempty(missing_cols)
+        @warn "Adding default entries for missing columns $(missing_cols)"
+    end
     # need to make sure columns are in the same order as in pc to avoid errors from TypedTables
     new_points = Table(NamedTuple{ (columnnames(pc)...,) }( (map(col -> hasproperty(points, col) ? getproperty(points, col) : zeros(eltype(getproperty(pc, col)), length(points)), columnnames(pc))...,) ))
     append!(pc, new_points)

From 0389032943ec91be30871da9829ed4db10d2c5c1 Mon Sep 17 00:00:00 2001
From: BenCurran98 <b.curran@fugro.com>
Date: Fri, 31 Jan 2025 13:33:34 +1000
Subject: [PATCH 16/16] Set legacy point count correctly to 0 when setting
 point format

Signed-off-by: BenCurran98 <b.curran@fugro.com>
---
 docs/src/interface.md |  8 +++++++-
 src/LASDatasets.jl    |  2 +-
 src/dataset.jl        | 17 +++++++++++++++++
 src/header.jl         | 14 ++++++++++++++
 test/dataset.jl       | 12 ++++++++++--
 5 files changed, 49 insertions(+), 4 deletions(-)

diff --git a/docs/src/interface.md b/docs/src/interface.md
index cb0ad7a..75afcae 100644
--- a/docs/src/interface.md
+++ b/docs/src/interface.md
@@ -77,12 +77,13 @@ save_las("my_las.las", pc)
 Note that when you supply just the point cloud outside of a `LASDataset`, *LASDatasets.jl* will automatically construct the appropriate header for you so you don't need to worry about the specifics of appropriate point formats etc. 
 
 ## Modifying LAS Contents
-You can modify point fields in your `LASDataset` by adding new columns or merging in values from an existing vector. Additionally, you can add new points to your dataset, however one should note that the user is responsible for correctly setting the appropriate fields such as synthetic flags themselves.
+You can modify point fields in your `LASDataset` by adding new columns or merging in values from an existing vector. Additionally, you can add new points to your dataset or remove them, however one should note when adding points that the user is responsible for correctly setting the appropriate fields such as synthetic flags themselves.
 
 ```@docs; canonical = false
 add_column!
 merge_column!
 add_points!
+remove_points!
 ```
 
 For example, if you want to add a set of synthetic points to your dataset, you can run:
@@ -94,6 +95,11 @@ synthetic_points = Table(position = rand(SVector{3, Float64}, 5), classification
 add_points!(las, synthetic_points)
 ```
 
+You can remove points from your data using the `remove_points!` function and specifying the indices of the points you wish to delete (these will be indexing into the list of points in order). E.g.
+```julia
+remove_points!(las, 11:15)
+```
+
 You can also add or remove *(E)VLRs* using the following functions, and set an existing *(E)VLR* as *superseded* if it's an old copy of a record.
 
 ```@docs; canonical = false
diff --git a/src/LASDatasets.jl b/src/LASDatasets.jl
index 6ba6f12..3a28805 100644
--- a/src/LASDatasets.jl
+++ b/src/LASDatasets.jl
@@ -55,7 +55,7 @@ export get_classes, get_description, set_description!
 export @register_vlr_type, read_vlr_data, extract_vlr_type
 
 export LASDataset, get_header, get_pointcloud, get_vlrs, get_evlrs, get_user_defined_bytes, get_unit_conversion
-export add_column!, add_points!, merge_column!, add_vlr!, remove_vlr!, set_superseded!
+export add_column!, add_points!, remove_points!, merge_column!, add_vlr!, remove_vlr!, set_superseded!
 
 # I/O methods
 export load_las, load_pointcloud, save_las, load_header, load_vlrs
diff --git a/src/dataset.jl b/src/dataset.jl
index abb59cf..bf9f097 100644
--- a/src/dataset.jl
+++ b/src/dataset.jl
@@ -455,6 +455,23 @@ function add_points!(las::LASDataset, points::AbstractVector{<:NamedTuple})
     return nothing
 end
 
+"""
+    $(TYPEDSIGNATURES)
+
+Remove a set of points stored at indices `idxs` from a `las` dataset. Updates header information to ensure consistency
+"""
+function remove_points!(las::LASDataset, idxs::Union{AbstractUnitRange, AbstractVector{<:Integer}})
+    pc = get_pointcloud(las)
+    @error "REMOVE"
+    @show length(pc)
+    @show length(pc.id)
+    deleteat!(pc, idxs)
+    @show length(pc)
+    @show length(pc.id)
+    _consolidate_point_header_info!(get_header(las), pc)
+    return nothing
+end
+
 # plumb through header util functions to act on a LASDataset for convenience
 for header_func ∈ (
     :las_version,
diff --git a/src/header.jl b/src/header.jl
index 694e616..5256f9d 100644
--- a/src/header.jl
+++ b/src/header.jl
@@ -593,6 +593,8 @@ function set_point_record_count!(header::LasHeader, num_points::Integer)
     header.record_count = UInt64(num_points)
     if get_point_format_id(point_format(header)) ≤ 5
         header.legacy_record_count = UInt32(num_points)
+    else
+        header.legacy_record_count = zero(UInt32)
     end
     return nothing
 end
@@ -616,6 +618,15 @@ function set_num_evlr!(header::LasHeader, n::Integer)
     header.n_evlr = UInt64(n)
 end
 
+"""
+    $(TYPEDSIGNATURES)
+
+Set the offset (in bytes) into the LAS file where the first EVLR occurs
+"""
+function set_evlr_start!(header, offset::Integer)
+    header.evlr_start = offset
+end
+
 """If true, GPS Time is standard GPS Time (satellite GPS Time) minus 1e9.
 If false, GPS Time is GPS Week Time.
 
@@ -873,7 +884,10 @@ end
 function _consolidate_point_header_info!(header::LasHeader, pointcloud::AbstractVector{<:NamedTuple})
     set_spatial_info!(header, get_spatial_info(pointcloud, scale(header)))
     
+    old_point_count = number_of_points(header)
     set_point_record_count!(header, length(pointcloud))
+    # make sure we update the EVLR start value if we change the point count
+    set_evlr_start!(header, evlr_start(header) + (length(pointcloud) - old_point_count))
     returns = (:returnnumber ∈ columnnames(pointcloud)) ? pointcloud.returnnumber : ones(Int, length(pointcloud))
     points_per_return = ntuple(r -> count(returns .== r), num_return_channels(header))
     set_number_of_points_by_return!(header, points_per_return)
diff --git a/test/dataset.jl b/test/dataset.jl
index 46e8abd..a5d444a 100644
--- a/test/dataset.jl
+++ b/test/dataset.jl
@@ -190,10 +190,9 @@
     add_column!(las, :overlap, falses(length(pc)))
     @test point_format(get_header(las)) == LasPoint6
     @test las_version(get_header(las)) == v"1.4"
-
     # we should be able to add new points in too
     new_points = Table(
-        id = length(pc):(length(pc) + 10),
+        id = (length(pc) + 1):(length(pc) + 10),
         position = rand(SVector{3, Float64}, 10),
         classification = rand(UInt8, 10),
         gps_time = rand(10),
@@ -207,8 +206,17 @@
     for col ∈ columnnames(new_points)
         @test getproperty(pointcloud, col)[101:end] == getproperty(new_points, col)
     end
+    orig_pc = deepcopy(pointcloud)
     # also make sure our header information is correctly set
     @test number_of_points(las) == 110
     @test spatial_info(las) == LASDatasets.get_spatial_info(pointcloud)
+
+    # we can also delete these points if we want
+    remove_points!(las, 101:110)
+    @test number_of_points(las) == 100
+    pc = get_pointcloud(las)
+    for col ∈ columnnames(pc)
+        @test getproperty(pc, col) == getproperty(orig_pc, col)[1:100]
+    end
     
 end
\ No newline at end of file