diff --git a/README.md b/README.md
index 792ac86..fc224db 100644
--- a/README.md
+++ b/README.md
@@ -20,7 +20,7 @@ table = Shapefile.Table(path)
 geoms = Shapefile.shapes(table)
 
 # whole columns can be retrieved by their name
-table.Descriptio  # => Union{String, Missing}["Square with triangle missing", "Smaller triangle", missing]
+table.Description  # => Union{String, Missing}["Square with triangle missing", "Smaller triangle", missing]
 
 # example function that iterates over the rows and gathers shapes that meet specific criteria
 function selectshapes(table)
@@ -50,6 +50,7 @@ julia> GeoInterface.coordinates(Shapefile.shape(first(table)))
 ```
 
 ## Alternative packages
+
 If you want another lightweight pure Julia package for reading feature files, consider
 also [GeoJSON.jl](https://github.com/JuliaGeo/GeoJSON.jl).
 
diff --git a/src/Shapefile.jl b/src/Shapefile.jl
index 06fea61..10981d6 100644
--- a/src/Shapefile.jl
+++ b/src/Shapefile.jl
@@ -12,6 +12,7 @@ include("polygons.jl")
 include("polylines.jl")
 include("multipoints.jl")
 include("multipatch.jl")
+include("utils.jl")
 
 # Handle/Table
 
diff --git a/src/common.jl b/src/common.jl
index 667cb23..dc8d98e 100644
--- a/src/common.jl
+++ b/src/common.jl
@@ -4,6 +4,11 @@ abstract type AbstractShape end
 GeoInterface.isgeometry(::Type{<:AbstractShape}) = true
 GeoInterface.ncoord(::GI.AbstractGeometryTrait, ::AbstractShape) = 2 # With specific methods when 3
 
+Base.convert(T::Type{<:AbstractShape}, geom) = Base.convert(T, GI.geomtrait(geom), geom)
+Base.convert(::Type{T}, geom::T) where {T<:AbstractShape} = geom
+Base.convert(::Type{<:AbstractShape}, ::Nothing, geom) =
+    throw(ArgumentError("object to convert is not a GeoInterface compatible geometry"))
+
 """
     Rect
 
diff --git a/src/extent.jl b/src/extent.jl
index 0afe46c..0552a02 100644
--- a/src/extent.jl
+++ b/src/extent.jl
@@ -1,14 +1,24 @@
-const Shape3D = Union{PolylineZ,PolygonZ,MultiPointZ,MultiPatch}
+const ShapeZ = Union{PolylineZ,PolygonZ,MultiPointZ,MultiPatch}
+const ShapeM = Union{ShapeZ,PolylineM,PolygonM,MultiPointM}
 
-# extent
+# 3d
 GI.is3d(::GI.AbstractGeometryTrait, ::AbstractShape) = false
-GI.is3d(::GI.AbstractGeometryTrait, ::Shape3D) = true
+GI.is3d(::GI.AbstractPointTrait, ::AbstractPoint) = false
+GI.is3d(::GI.AbstractGeometryTrait, ::ShapeZ) = true
+GI.is3d(::GI.AbstractPointTrait, ::PointZ) = true
+
+# measured
+GI.ismeasured(::GI.AbstractGeometryTrait, ::AbstractShape) = false
+GI.ismeasured(::GI.AbstractGeometryTrait, ::ShapeM) = true
+GI.ismeasured(::GI.AbstractPointTrait, ::Point) = false
+GI.ismeasured(::GI.AbstractPointTrait, ::Union{PointM,PointZ}) = true
 
+# extent
 function GI.extent(x::AbstractShape)
     rect = x.MBR
     return Extents.Extent(X=(rect.left, rect.right), Y=(rect.bottom, rect.top))
 end
-function GI.extent(x::Union{Shape3D,Handle})
+function GI.extent(x::Union{ShapeZ,Handle})
     rect = x.MBR
     return Extents.Extent(X=(rect.left, rect.right), Y=(rect.bottom, rect.top), Z=(x.zrange.left, x.zrange.right))
 end
diff --git a/src/multipatch.jl b/src/multipatch.jl
index 2795dcb..d125b48 100644
--- a/src/multipatch.jl
+++ b/src/multipatch.jl
@@ -40,15 +40,11 @@ function Base.read(io::IO, ::Type{MultiPatch})
     box = read(io, Rect)
     numparts = read(io, Int32)
     numpoints = read(io, Int32)
-    parts = Vector{Int32}(undef, numparts)
-    read!(io, parts)
+    parts = _readparts(io, numparts)
     parttypes = Vector{Int32}(undef, numparts)
     read!(io, parttypes)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    zrange = read(io, Interval)
-    zvalues = Vector{Float64}(undef, numpoints)
-    read!(io, zvalues)
+    points = _readpoints(io, numpoints)
+    zrange, zvalues = _readfloats(io, numpoints)
     # mrange = Vector{Float64}(2)
     # read!(io, mrange)
     # measures = Vector{Float64}(numpoints)
diff --git a/src/multipoints.jl b/src/multipoints.jl
index 32e89df..80d5542 100644
--- a/src/multipoints.jl
+++ b/src/multipoints.jl
@@ -1,6 +1,11 @@
 
 abstract type AbstractMultiPoint{T} <: AbstractShape end
 
+_pointtype(::Type{<:AbstractMultiPoint{T}}) where T = T
+
+Base.convert(::Type{T}, ::GI.MultiPointTrait, geom) where T<:AbstractMultiPoint =
+    T(_convert(_pointtype(T), geom)...)
+
 GI.geomtrait(::AbstractMultiPoint) = GI.MultiPointTrait()
 GI.ngeom(::GI.MultiPointTrait, geom::AbstractMultiPoint) = length(geom.points)
 GI.ncoord(::GI.MultiPointTrait, geom::AbstractMultiPoint{T}) where {T} = _ncoord(T)
@@ -23,13 +28,13 @@ end
 function Base.read(io::IO, ::Type{MultiPoint})
     box = read(io, Rect)
     numpoints = read(io, Int32)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    MultiPoint(box, points)
+    points = _readpoints(io, numpoints)
+    return MultiPoint(box, points)
 end
 
 GI.getgeom(::GI.MultiPointTrait, geom::MultiPoint, i::Integer) = geom.points[i]
 
+
 """
     MultiPointM <: AbstractMultiPoint
 
@@ -54,11 +59,8 @@ end
 function Base.read(io::IO, ::Type{MultiPointM})
     box = read(io, Rect)
     numpoints = read(io, Int32)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    mrange = read(io, Interval)
-    measures = Vector{Float64}(undef, numpoints)
-    read!(io, measures)
+    points = _readpoints(io, numpoints)
+    mrange, measures = _readm(io, numpoints)
     MultiPointM(box, points, mrange, measures)
 end
 
@@ -95,13 +97,8 @@ GI.getgeom(::GI.MultiPointTrait, geom::MultiPointZ, i::Integer) =
 function Base.read(io::IO, ::Type{MultiPointZ})
     box = read(io, Rect)
     numpoints = read(io, Int32)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    zrange = read(io, Interval)
-    zvalues = Vector{Float64}(undef, numpoints)
-    read!(io, zvalues)
-    mrange = read(io, Interval)
-    measures = Vector{Float64}(undef, numpoints)
-    read!(io, measures)
+    points = _readpoints(io, numpoints)
+    zrange, zvalues = _readz(io, numpoints)
+    mrange, measures = _readm(io, numpoints)
     MultiPointZ(box, points, zrange, zvalues, mrange, measures)
 end
diff --git a/src/points.jl b/src/points.jl
index fb324ab..a7a6786 100644
--- a/src/points.jl
+++ b/src/points.jl
@@ -27,6 +27,8 @@ function Base.read(io::IO, ::Type{Point})
     Point(x, y)
 end
 
+Base.convert(::Type{<:Point}, ::GI.PointTrait, geom) = Point(GI.x(geom), GI.y(geom))
+
 """
     PointM <: AbstractPoint
 
@@ -50,6 +52,11 @@ function Base.read(io::IO, ::Type{PointM})
     PointM(x, y, m)
 end
 
+function Base.convert(::Type{<:PointM}, ::GI.PointTrait, geom)
+    m = GI.ismeasured(geom) ? GI.m(geom) : 0.0
+    PointM(GI.x(geom), GI.y(geom), m)
+end
+
 GI.m(::GI.PointTrait, point::PointM) = point.m
 
 """
@@ -77,6 +84,12 @@ function Base.read(io::IO, ::Type{PointZ})
     PointZ(x, y, z, m)
 end
 
+function Base.convert(::Type{<:PointZ}, ::GI.PointTrait, geom)
+    z = GI.is3d(geom) ? GI.z(geom) : 0.0
+    m = GI.ismeasured(geom) ? GI.m(geom) : 0.0
+    PointZ(GI.x(geom), GI.y(geom), z, m)
+end
+
 GI.m(::GI.PointTrait, point::PointZ) = point.m
 GI.z(::GI.PointTrait, point::PointZ) = point.z
 
diff --git a/src/polygons.jl b/src/polygons.jl
index 8e6d5bd..ca41cc1 100644
--- a/src/polygons.jl
+++ b/src/polygons.jl
@@ -23,7 +23,7 @@ GI.getgeom(::GI.LinearRingTrait, lr::LinearRing, i::Integer) = lr[i]
 
 Base.getindex(lr::LinearRing{Point}, i) = lr.xy[i]
 Base.getindex(lr::LinearRing{PointM}, i) = PointM(lr.xy[i], lr.m[i])
-Base.getindex(lr::LinearRing{PointZ}, i) = PointZ(lr.xy[i], lr.m[i], lr.z[i])
+Base.getindex(lr::LinearRing{PointZ}, i) = PointZ(lr.xy[i], lr.z[i], lr.m[i])
 Base.size(lr::LinearRing) = (length(lr),)
 Base.length(lr::LinearRing) = length(lr.xy)
 
@@ -45,6 +45,11 @@ Base.push!(p::SubPolygon, x) = Base.push!(parent(p), x)
 
 abstract type AbstractPolygon{T} <: AbstractShape end
 
+_pointtype(::Type{<:AbstractPolygon{T}}) where T = T
+
+Base.convert(::Type{T}, ::GI.MultiPolygonTrait, geom) where T<:AbstractPolygon =
+    T(_convertparts(_pointtype(T), geom)...)
+
 # Shapefile polygons are OGC multipolygons
 GI.geomtrait(geom::AbstractPolygon) = GI.MultiPolygonTrait()
 GI.nring(::GI.MultiPolygonTrait, geom::AbstractPolygon) = length(geom.parts)
@@ -59,6 +64,7 @@ function GI.ngeom(::GI.MultiPolygonTrait, geom::AbstractPolygon)
     return n
 end
 
+
 function GI.getring(t::GI.MultiPolygonTrait, geom::AbstractPolygon)
     return (GI.getring(t, geom, i) for i in eachindex(geom.parts))
 end
@@ -146,10 +152,8 @@ function Base.read(io::IO, ::Type{Polygon})
     box = read(io, Rect)
     numparts = read(io, Int32)
     numpoints = read(io, Int32)
-    parts = Vector{Int32}(undef, numparts)
-    read!(io, parts)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
+    parts = _readparts(io, numparts)
+    points = _readpoints(io, numpoints)
     Polygon(box, parts, points)
 end
 
@@ -176,13 +180,9 @@ function Base.read(io::IO, ::Type{PolygonM})
     box = read(io, Rect)
     numparts = read(io, Int32)
     numpoints = read(io, Int32)
-    parts = Vector{Int32}(undef, numparts)
-    read!(io, parts)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    mrange = read(io, Interval)
-    measures = Vector{Float64}(undef, numpoints)
-    read!(io, measures)
+    parts = _readparts(io, numparts)
+    points = _readpoints(io, numpoints)
+    mrange, measures = _readm(io, numpoints)
     PolygonM(box, parts, points, mrange, measures)
 end
 
@@ -212,15 +212,9 @@ function Base.read(io::IO, ::Type{PolygonZ})
     box = read(io, Rect)
     numparts = read(io, Int32)
     numpoints = read(io, Int32)
-    parts = Vector{Int32}(undef, numparts)
-    read!(io, parts)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    zrange = read(io, Interval)
-    zvalues = Vector{Float64}(undef, numpoints)
-    read!(io, zvalues)
-    mrange = read(io, Interval)
-    measures = Vector{Float64}(undef, numpoints)
-    read!(io, measures)
+    parts = _readparts(io, numparts)
+    points = _readpoints(io, numpoints)
+    zrange, zvalues = _readz(io, numpoints)
+    mrange, measures = _readm(io, numpoints)
     PolygonZ(box, parts, points, zrange, zvalues, mrange, measures)
 end
diff --git a/src/polylines.jl b/src/polylines.jl
index 493cd97..ced4945 100644
--- a/src/polylines.jl
+++ b/src/polylines.jl
@@ -19,11 +19,16 @@ Base.size(p::LineString) = (length(p),)
 Base.length(lr::LineString) = length(parent(lr))
 Base.getindex(lr::LineString{Point}, i) = lr.xy[i]
 Base.getindex(lr::LineString{PointM}, i) = PointM(lr.xy[i], lr.m[i])
-Base.getindex(lr::LineString{PointZ}, i) = PointZ(lr.xy[i], lr.m[i], lr.z[i])
+Base.getindex(lr::LineString{PointZ}, i) = PointZ(lr.xy[i], lr.z[i], lr.m[i])
 
 
 abstract type AbstractPolyline{T} <: AbstractShape end
 
+_pointtype(::Type{<:AbstractPolyline{T}}) where T = T
+
+Base.convert(::Type{T}, ::GI.MultiLineStringTrait, geom) where T<:AbstractPolyline =
+    T(_convertparts(_pointtype(T), geom)...)
+
 GI.geomtrait(::AbstractPolyline) = GI.MultiLineStringTrait()
 GI.ngeom(::GI.MultiLineStringTrait, geom::AbstractPolyline) = length(geom.parts)
 GI.ncoord(::GI.MultiLineStringTrait, ::AbstractPolyline{T}) where {T} = _ncoord(T)
@@ -59,10 +64,8 @@ function Base.read(io::IO, ::Type{Polyline})
     box = read(io, Rect)
     numparts = read(io, Int32)
     numpoints = read(io, Int32)
-    parts = Vector{Int32}(undef, numparts)
-    read!(io, parts)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
+    parts = _readparts(io, numparts)
+    points = _readpoints(io, numpoints)
     Polyline(box, parts, points)
 end
 
@@ -94,13 +97,9 @@ function Base.read(io::IO, ::Type{PolylineM})
     box = read(io, Rect)
     numparts = read(io, Int32)
     numpoints = read(io, Int32)
-    parts = Vector{Int32}(undef, numparts)
-    read!(io, parts)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    mrange = read(io, Interval)
-    measures = Vector{Float64}(undef, numpoints)
-    read!(io, measures)
+    parts = _readparts(io, numparts)
+    points = _readpoints(io, numpoints)
+    mrange, measures = _readm(io, numpoints)
     PolylineM(box, parts, points, mrange, measures)
 end
 
@@ -133,15 +132,9 @@ function Base.read(io::IO, ::Type{PolylineZ})
     box = read(io, Rect)
     numparts = read(io, Int32)
     numpoints = read(io, Int32)
-    parts = Vector{Int32}(undef, numparts)
-    read!(io, parts)
-    points = Vector{Point}(undef, numpoints)
-    read!(io, points)
-    zrange = read(io, Interval)
-    zvalues = Vector{Float64}(undef, numpoints)
-    read!(io, zvalues)
-    mrange = read(io, Interval)
-    measures = Vector{Float64}(undef, numpoints)
-    read!(io, measures)
+    parts = _readparts(io, numparts)
+    points = _readpoints(io, numpoints)
+    zrange, zvalues = _readz(io, numpoints)
+    mrange, measures = _readm(io, numpoints)
     PolylineZ(box, parts, points, zrange, zvalues, mrange, measures)
 end
diff --git a/src/utils.jl b/src/utils.jl
new file mode 100644
index 0000000..8948acc
--- /dev/null
+++ b/src/utils.jl
@@ -0,0 +1,100 @@
+_partvec(n) = Vector{Int32}(undef, n)
+_pointvec(n) = Vector{Point}(undef, n)
+_floatvec(n) = Vector{Float64}(undef, n)
+
+function _readparts(io, n)
+    points = _partvec(n)
+    read!(io, points)
+    return points
+end
+
+function _readpoints(io, n)
+    points = _pointvec(n)
+    read!(io, points)
+    return points
+end
+
+_readm(io, n) = _readfloats(io, n)
+_readz(io, n) = _readfloats(io, n)
+function _readfloats(io, n)
+    interval = read(io, Interval)
+    values = _floatvec(n)
+    read!(io, values)
+    return interval, values
+end
+
+# Most objects have similar structure, so we share conversion internals
+function _convert(pointmode::Type{<:Point}, geom)
+    n = GI.npoint(geom)
+    points = _pointvec(n)
+    for (i, point) in enumerate(GI.getpoint(geom))
+        points[i] = Point(GI.x(point), GI.y(point))
+    end
+    MBR = _getbounds(points)
+    return MBR, points
+end
+function _convert(::Type{<:PointM}, geom)
+    hasm = GI.ismeasured(first(GI.getpoint(geom)))
+    n = GI.npoint(geom)
+    points = _pointvec(n)
+    measures = _floatvec(n)
+    for (i, point) in enumerate(GI.getpoint(geom))
+        points[i] = Point(GI.x(point), GI.y(point))
+        measures[i] = hasm ? GI.m(point) : 0.0
+    end
+    @show measures geom.measures
+    mrange = Interval(extrema(measures)...)
+    MBR = _getbounds(points)
+    return MBR, points, mrange, measures
+end
+function _convert(::Type{<:PointZ}, geom)
+    hasz = GI.is3d(first(GI.getpoint(geom))) 
+    hasm = GI.ismeasured(first(GI.getpoint(geom)))
+    n = GI.npoint(geom)
+    points = _pointvec(n)
+    zvalues = _floatvec(n)
+    measures = _floatvec(n)
+    for (i, point) in enumerate(GI.getpoint(geom))
+        points[i] = Point(GI.x(point), GI.y(point))
+        zvalues[i] = hasz ? GI.z(point) : 0.0
+        measures[i] = hasm ? GI.m(point) : 0.0
+    end
+    mrange = Interval(extrema(measures)...)
+    zrange = Interval(extrema(zvalues)...)
+    MBR = _getbounds(points)
+    return MBR, points, zrange, zvalues, mrange, measures
+end
+
+# _convert with additional parts field
+function _convertparts(T, geom) 
+    MBR, args... = _convert(T, geom)
+    parts = _getparts(geom)
+    return (MBR, parts, args...)
+end
+
+_getparts(geom) = _getparts(GI.geomtrait(geom), geom)
+# Special-case multi polygons because we only store the rings
+function _getparts(::GI.MultiPolygonTrait, geom)
+    parts = _partvec(GI.nring(geom))
+    n = 0
+    for (i, ring) in enumerate(GI.getring(geom))
+        parts[i] = n
+        n += GI.npoint(ring)
+    end
+    return parts
+end
+function _getparts(::GI.AbstractGeometryTrait, geom)
+    parts = _partvec(GI.ngeom(geom))
+    n = 0
+    for (i, g) in enumerate(GI.getgeom(geom))
+        parts[i] = n
+        n += GI.npoint(g)
+    end
+    return parts
+end
+
+function _getbounds(points::Vector{Point})
+    xrange = extrema(GI.x(p) for p in points)
+    yrange = extrema(GI.y(p) for p in points)
+    return Rect(xrange[1], yrange[1], xrange[2], yrange[2])
+end
diff --git a/test/runtests.jl b/test/runtests.jl
index 5a4b830..37d11f6 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -98,29 +98,61 @@ test_tuples = [
 #     end
 # end
 
+points = [Point(1, 3), Point(2, 3), Point(2, 4), Point(1, 4)]
+test_shapes = Dict(
+    Point => Point(1, 2),
+    PointM => PointM(1, 2, 3),
+    PointZ => PointZ(1, 2, 3, 4),
+    MultiPoint => MultiPoint(Rect(1, 3, 2, 4), points),
+    MultiPointM => MultiPointM(Rect(1, 3, 2, 4), points, Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
+    MultiPointZ => MultiPointZ(Rect(1, 3, 2, 4), points, Interval(1, 4), [1, 2, 3, 4], Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
+    Polygon => Polygon(Rect(1, 3, 2, 4), [0], points),
+    PolygonM => PolygonM(Rect(1, 3, 2, 4), [0], points, Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
+    PolygonZ => PolygonZ(Rect(1, 3, 2, 4), [0], points, Interval(1, 4), [1, 2, 3, 4], Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
+    Polyline => Polyline(Rect(1, 3, 2, 4), [0], points),
+    PolylineM => PolylineM(Rect(1, 3, 2, 4), [0], points, Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
+    PolylineZ => PolylineZ(Rect(1, 3, 2, 4), [0], points, Interval(1, 4), [1, 2, 3, 4], Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
+    LineString => LineString{Point}(view(points, 1:4)),
+    LinearRing => LinearRing{Point}(view(points, 1:4)),
+    SubPolygon => SubPolygon([LinearRing{Point}(view(points, 1:4))]),
+)
+
 @testset "GeoInterface compatability" begin
-    points = [Point(1, 2), Point(2, 2), Point(2, 1), Point(1, 1)]
-
-    shapes = (
-        Point(1, 2),
-        PointM(1, 2, 3),
-        PointZ(1, 2, 3, 4),
-        MultiPoint(Rect(1, 2, 2, 1), points),
-        MultiPointM(Rect(1, 2, 2, 1), points, Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
-        MultiPointZ(Rect(1, 2, 2, 1), points, Interval(1, 4), [1, 2, 3, 4], Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
-        Polygon(Rect(1, 2, 2, 1), [1], points),
-        PolygonM(Rect(1, 2, 2, 1), [1], points, Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
-        PolygonZ(Rect(1, 2, 2, 1), [1], points, Interval(1, 4), [1, 2, 3, 4], Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
-        Polyline(Rect(1, 2, 2, 1), [1], points),
-        PolylineM(Rect(1, 2, 2, 1), [1], points, Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
-        PolylineZ(Rect(1, 2, 2, 1), [1], points, Interval(1, 4), [1, 2, 3, 4], Interval(10.0, 40.0), [10.0, 20.0, 30.0, 40.0]),
-        LineString{Point}(view(points, 1:4)),
-        LinearRing{Point}(view(points, 1:4)),
-        SubPolygon([LinearRing{Point}(view(points, 1:4))]),
-    )
+    @test all(s -> GeoInterface.testgeometry(s), values(test_shapes))
+    @test_broken GeoInterface.testgeometry(MultiPatch(Rect(1, 3, 2, 4), [0], [1], points, Interval(1, 4), [1, 2, 3, 4]))
+end
 
-    @test all(s -> GeoInterface.testgeometry(s), shapes)
-    @test_broken GeoInterface.testgeometry(MultiPatch(Rect(1, 2, 2, 1), [1], [1], points, Interval(1, 4), [1, 2, 3, 4]))
+@testset "convert" begin
+    @test convert(Point, PointZ(1, 2, 3, 4)) == Point(1, 2)
+    @test convert(PointM, PointZ(1, 2, 3, 4)) == PointM(1, 2, 4)
+    @test convert(PointZ, Point(1, 2)) == PointZ(1, 2, 0, 0)
+    pl = convert(Polyline, test_shapes[PolylineZ])
+    @test pl.MBR == test_shapes[Polyline].MBR
+    @test pl.points == test_shapes[Polyline].points
+    pl = convert(PolylineM, test_shapes[PolylineZ])
+    @test pl.MBR == test_shapes[PolylineM].MBR
+    @test pl.points == test_shapes[PolylineM].points
+    @test pl.measures == test_shapes[PolylineM].measures
+    pg = convert(Polygon, test_shapes[PolygonZ])
+    @test pg.MBR == test_shapes[Polygon].MBR
+    @test pg.points == test_shapes[Polygon].points
+    pg = convert(PolygonM, test_shapes[PolygonZ])
+    @test pg.MBR == test_shapes[PolygonM].MBR
+    @test pg.points == test_shapes[PolygonM].points
+    @test pg.measures == test_shapes[PolygonM].measures
+    mp = convert(MultiPoint, test_shapes[MultiPointZ])
+    @test mp.MBR == test_shapes[MultiPoint].MBR
+    @test mp.points == test_shapes[MultiPoint].points
+    mp = convert(MultiPointM, test_shapes[MultiPointZ])
+    @test mp.MBR == test_shapes[PolylineM].MBR
+    @test mp.points == test_shapes[PolylineM].points
+    @test mp.measures == test_shapes[PolylineM].measures
+
+    mp = convert(MultiPointZ, test_shapes[MultiPointM])
+    @test mp.MBR == test_shapes[PolylineZ].MBR
+    @test mp.points == test_shapes[PolylineZ].points
+    @test mp.measures == test_shapes[PolylineM].measures
+    @test mp.zvalues == [0.0, 0.0, 0.0, 0.0]  # Missing Z filled with zeros: is this the best default?
 end
 
 @testset "Loading Shapefiles" begin