use StepRange rather than LinRange

ext/RastersArchGDALExt/gdal_source.jl

@@ -174,8 +174,7 @@ function DD.dims(raster::AG.RasterDataset, crs=nothing, mappedcrs=nothing)
             xmax = gt[GDAL_TOPLEFT_X] + xstep * xsize
             xorder = ReverseOrdered()
         end
-        xindex = LinRange(xmin, xmax, xsize)
-        xorder = xstep > 0 ? ForwardOrdered() : ReverseOrdered()
+        xindex = xmin:xstep:xmax 
 
         ystep = gt[GDAL_NS_RES] # A negative number
         if ystep > 0 
@@ -187,7 +186,7 @@ function DD.dims(raster::AG.RasterDataset, crs=nothing, mappedcrs=nothing)
             ymin = gt[GDAL_TOPLEFT_Y] + ystep * ysize
             yorder = ReverseOrdered()
         end
-        yindex = LinRange(ymax, ymin, ysize)
+        yindex = ymax:ystep:ymin
 
         # Spatial data defaults to area/inteval
         xsampling, ysampling = if _gdalmetadata(raster.ds, "AREA_OR_POINT") == "Point"
@@ -199,7 +198,7 @@ function DD.dims(raster::AG.RasterDataset, crs=nothing, mappedcrs=nothing)
 
         xlookup = Projected(xindex;
             order=xorder,
-            span=Regular(step(xindex)),
+            span=Regular(xstep),
             sampling=xsampling,
             metadata=xy_metadata,
             crs=crs,
@@ -209,7 +208,7 @@ function DD.dims(raster::AG.RasterDataset, crs=nothing, mappedcrs=nothing)
             order=yorder,
             sampling=ysampling,
             # Use the range step as is will be different to ystep due to float error
-            span=Regular(step(yindex)),
+            span=Regular(ystep),
             metadata=xy_metadata,
             crs=crs,
             mappedcrs=mappedcrs,

src/lookup.jl

@@ -184,13 +184,13 @@ end
 
 
 
-_projectedrange(l::Projected) = LinRange(first(l), last(l), length(l))
+_projectedrange(l::Projected) = first(l):step(span(l)):last(l)
 _projectedrange(l::Mapped) = _projectedrange(span(l), crs(l), l)
 function _projectedrange(span, crs, l::Mapped)
     start, stop = reproject(mappedcrs(l), crs, dim(l), [first(l), last(l)])
-    LinRange(start, stop, length(l))
+    start:(stop - start) / length(l):stop
 end
-_projectedrange(::Regular, crs::Nothing, l::Mapped) = LinRange(first(l), last(l), length(l))
+_projectedrange(::Regular, crs::Nothing, l::Mapped) = first(l):step(span(l)):last(l)
 function _projectedrange(::T, crs::Nothing, l::Mapped) where T<:Union{Irregular,Explicit}
     error("Cannot convert a Mapped $T index to Projected when crs is nothing")
 end

src/methods/aggregate.jl

@@ -232,7 +232,7 @@ function disaggregate(locus, lookup::LookupArray, scale)
     step_ = step(lookup) / intscale
     start = lookup[1] - _agoffset(locus, intscale) * step_
     stop = start + (len - 1)  * step_
-    index = LinRange(start, stop, len)
+    index = start:step_:stop
     if lookup isa AbstractSampled
         sp = disaggregate(locus, span(lookup), intscale)
         return rebuild(lookup; data=index, span=sp)

src/sources/grd.jl

@@ -50,18 +50,18 @@ function DD.dims(grd::GRDattrib, crs=nothing, mappedcrs=nothing)
     ybounds = parse.(Float64, (attrib["ymin"], attrib["ymax"]))
 
     # Always intervals
-    xspan = (xbounds[2] - xbounds[1]) / xsize
-    yspan = (ybounds[1] - ybounds[2]) / ysize
+    xstep = (xbounds[2] - xbounds[1]) / xsize
+    ystep = (ybounds[1] - ybounds[2]) / ysize
 
     # Not fully implemented yet
     xy_metadata = _metadatadict(GRDsource)
 
-    xindex = LinRange(xbounds[1], xbounds[2] - xspan, xsize)
-    yindex = LinRange(ybounds[2] + yspan, ybounds[1], ysize)
+    xindex = xbounds[1] : xstep : xbounds[2] - xstep
+    yindex = ybounds[2] + ystep : ystep : ybounds[1]
 
     xlookup = Projected(xindex;
         order=GRD_X_ORDER,
-        span=Regular(xspan),
+        span=Regular(xstep),
         sampling=Intervals(Start()),
         metadata=xy_metadata,
         crs=crs,
@@ -70,7 +70,7 @@ function DD.dims(grd::GRDattrib, crs=nothing, mappedcrs=nothing)
     )
     ylookup = Projected(yindex;
         order=GRD_Y_ORDER,
-        span=Regular(yspan),
+        span=Regular(ystep),
         sampling=Intervals(Start()),
         metadata=xy_metadata,
         crs=crs,

test/aggregate.jl

@@ -26,8 +26,8 @@ data1 = [ 1  2  3  4  5  6 -1
 data2 = 2 * data1
 data3 = 3 * data1
 data4 = 4 * data1
-dimz = X(Sampled([30., 40., 50.]; order=ForwardOrdered(), span=Regular(10.0), sampling=Points())), 
-       Y(Sampled(LinRange(-10., 20., 7); order=ForwardOrdered(), span=Regular(5.0), sampling=Points()))
+dimz = X(Sampled([30.0, 40.0, 50.]; order=ForwardOrdered(), span=Regular(10.0), sampling=Points())), 
+Y(Sampled(-10.0:30//6:20.0; order=ForwardOrdered(), span=Regular(5.0), sampling=Points()))
 array1 = Raster(data1, dimz)
 array2 = Raster(data2, dimz)
 array1a = Raster(data3, dimz)
@@ -39,7 +39,7 @@ series = RasterSeries([stack1, stack2], (Ti(dates),));
 
 
 @testset "Aggregate a dimension" begin
-    lat = Y(Sampled(LinRange(3, 13, 6), ForwardOrdered(), Regular(2.0), Intervals(Start()), NoMetadata()))
+    lat = Y(Sampled(3.0:2.0:13.0, ForwardOrdered(), Regular(2.0), Intervals(Start()), NoMetadata()))
     aglat = aggregate(Start(), lat, 3)
     @test span(lookup(aglat)) == Regular(6.0)
     @test disaggregate(Start(), aglat, 3) == lat

test/reproject.jl

@@ -37,7 +37,7 @@ end
 
     lonstart, lonend = 0.5, 179.5
     cealonstart, cealonend = reproject(proj4326, projcea, X(), [lonstart, lonend])
-    cealonrange = LinRange(cealonstart, cealonend, 180)
+    cealonrange = cealonstart:1.0:cealonend
     lon = X(Projected(cealonrange, ForwardOrdered(), Regular(step(cealonrange)), 
                       Intervals(Center()), NoMetadata(), projcea, proj4326, X()))
     convertedlon = convertlookup(Mapped, lon)
@@ -52,7 +52,7 @@ end
     latbounds = -90.0, 90.0
     ceabounds = reproject(proj4326, projcea, Y(), latbounds)
     # Vals are the bounding range without the end bound
-    cealatrange = LinRange(ceabounds..., 181)[1:180]
+    cealatrange = ceabounds[1] : (ceabounds[2] - ceabounds[1] / 180) : ceabounds[2]
     lat = Y(Projected(cealatrange, ForwardOrdered(), Regular(step(cealatrange)), 
                       Intervals(Start()), NoMetadata(), projcea, proj4326, Y()))
     convertedlat = convertlookup(Mapped, lat)

test/sources/gdal.jl

@@ -479,7 +479,7 @@ gdalpath = maybedownload(url)
         @test order(dims(rast)) == (ForwardOrdered(), ForwardOrdered())
         @test span(rast) == (Regular(1.0), Regular(1.0))
         @test sampling(rast) == (Intervals(Start()), Intervals(Start()))
-        @test index(rast) == (LinRange(0.0, 239.0, 240), LinRange(0.0, 179.0, 180))
+        @test index(rast) == (0.0::1.0:239.0, 0.0:1.0:179.0)
     end
 
 end