updates

src/methods/crop_extend.jl

@@ -172,27 +172,25 @@ function _extend_to(x::RasterStackOrArray, to; kw...)
 end
 
 function _extend_to(A::AbstractRaster, to::DimTuple;
-    filename=nothing, suffix=nothing, touches=false
+    filename=nothing, suffix=nothing, touches=false, missingval=missingval(A)
 )
     others = otherdims(to, A)
     # Allow not specifying all dimensions
     to = (set(dims(A), map(=>, dims(A, to), to)...)..., others...)
     # Calculate the range of the old array in the extended array
-    ranges = _without_mapped_crs(A) do A
+    rangedims = _without_mapped_crs(A) do A
         _without_mapped_crs(to) do to
-            map(dims(A), to) do d, t
-                range = if touches 
-                    DD.selectindices(t, LA.Touches(bounds(d)))
-                else
-                    DD.selectindices(t, LA.ClosedInterval(bounds(d)...))
-                end
-                rebuild(d, range)
+            map(dims(A, to), to) do d, t
+                # Values must match exactly, so use `At`
+                DD.selectindices(t, At(first(d))):DD.selectindices(t, At(last(d)))
             end
         end
     end
+    others = otherdims(to, A)
+    final_to = (set(dims(A), map(=>, dims(A, to), to)...)..., others...)
     # Create a new extended array
-    newA = create(filename, eltype(A), to;
-        suffix, parent=parent(A), missingval=missingval(A),
+    newA = create(filename, eltype(A), final_to;
+        suffix, parent=parent(A), missingval,
         name=name(A), metadata=metadata(A)
     )
     # Input checks
@@ -208,15 +206,15 @@ function _extend_to(A::AbstractRaster, to::DimTuple;
         end
     end
     # The missingval may have changed for disk-based arrays
-    if !isequal(missingval(A), missingval(newA))
-        A = replace_missing(A, missingval(newA))
+    if !isequal(missingval, Rasters.missingval(newA))
+        A = replace_missing(A, Rasters.missingval(newA))
     end
     open(newA; write=true) do O
         # Fill it with missing/nodata values
-        O .= missingval(O)
+        O .= Rasters.missingval(O)
         # Copy the original data to the new array
         # Somehow this is slow from disk?
-        broadcast_dims!(identity, view(O, ranges...), A)
+        broadcast_dims!(identity, view(O, rangedims...), A)
     end
     return newA
 end