add test for merge_var

src/sf/bbox.jl

@@ -24,6 +24,12 @@ function bbox2cellsize(b::bbox, size)
 end
 
 # 两个至少提供一个
+"""
+    bbox2dims(b::bbox; size=nothing, cellsize=nothing, reverse_lat=true)
+
+# Arguments
+- `cellsize`: celly为负，则lat倒序
+"""
 function bbox2dims(b::bbox; size=nothing, cellsize=nothing, reverse_lat=true)
   if size !== nothing && cellsize === nothing
     cellsize = bbox2cellsize(b, size)

src/sf/st_dims.jl

@@ -48,6 +48,8 @@ function st_cellsize(lon::AbstractVector, lat::AbstractVector)
   lon[2] - lon[1], lat[2] - lat[1] # cellx, celly
 end
 
+st_cellsize(f::AbstractString) = gdalinfo(f)["cellsize"]
+
 # function st_cellsize(ra::AbstractRaster)
 #   # x, y = st_dims(r)
 #   # x[2] - x[1], y[2] - y[1]

src/sf/st_mosaic.jl

@@ -32,32 +32,40 @@ function st_mosaic(rs::Vector{SpatRaster{T,N}}; kw...) where {T,N}
 end
 
 
-function merge_var!(R::AbstractArray, f; var=nothing, box::bbox)
+function merge_var!(R::AbstractArray, f;
+  var=nothing, box::bbox, cellsize=nothing)
+
+  isnothing(cellsize) && (cellsize = st_cellsize(f))
   b = st_bbox(f)
   bands = bandnames(f)
   nband = length(bands)
   inds_var = !isnothing(var) ? grep(bands, var) : 1:nband
   ntime = length(inds_var)
 
   println("Reading data ...")
-  @time A = read_gdal(f, inds_var)
+  A = read_gdal(f, inds_var)
   ilon, ilat = bbox_overlap(b, box; cellsize)
   R[ilon, ilat, 1:ntime] .= A
+  nothing
 end
 
 function merge_var(fs; vars=nothing, var=nothing,
+  progress=true,
   box::bbox=bbox(-180, -60, 180, 90))
 
   f = fs[1]
-  cellsize = gdalinfo(f)["cellsize"][1]
+  cellsize = st_cellsize(f)
   lon, lat = bbox2dims(box; cellsize)
   nlon, nlat = length(lon), length(lat)
+
   bands = bandnames(f)
   ntime = isnothing(vars) ? length(bands) : length(bands) / length(vars)
 
   R = zeros(Float32, nlon, nlat, ntime)
-  @showprogress for f in fs
-    merge_var!(R, f; var, box)
+  p = Progress(ntime)
+  for f in fs
+    progress && next!(p)
+    merge_var!(R, f; var, box, cellsize)
   end
   R
 end

test/sf/test-st_mosaic.jl

@@ -1,10 +1,21 @@
-@testset "st_mosaic" begin  
+using Test, Ipaper, Ipaper.sf, ArchGDAL
+
+# box = st_bbox(st_bbox.(fs))
+@testset "st_mosaic" begin
   bands = string.(1:4)
-  r1 = rast(rand(4, 4, 4), bbox(-180.0, -60.0, 180.0, -30.0); bands)
   r2 = rast(rand(4, 4, 4), bbox(-180.0, -30.0, 180.0, 0.0); bands)
-
+  r1 = rast(rand(4, 4, 4), bbox(-180.0, -60.0, 180.0, -30.0); bands)
+
   rs = [r1, r2]
   r_big = st_mosaic(rs)
-  @test st_bbox(r_big) == bbox(-180.0, -60.0, 180.0, 0.0) 
-  # r_big
+  @test st_bbox(r_big) == bbox(-180.0, -60.0, 180.0, 0.0)
+
+  ## test for `merge_var`
+  fs = ["r1.tif", "r2.tif"]
+  write_gdal(r1, fs[1])
+  write_gdal(r2, fs[2])
+
+  R = merge_var(fs; box)
+  @test cat(r2.A, r1.A, dims=2) ≈ R
+  rm.(fs)
 end