diff --git a/Project.toml b/Project.toml index 26d1b470..083fa7f9 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "SimpleSDMLayers" uuid = "2c645270-77db-11e9-22c3-0f302a89c64c" authors = ["Timothée Poisot ", "Gabriel Dansereau "] -version = "0.4.8" +version = "0.4.9" [deps] ArchGDAL = "c9ce4bd3-c3d5-55b8-8973-c0e20141b8c3" diff --git a/src/SimpleSDMLayers.jl b/src/SimpleSDMLayers.jl index 0c2d1ae2..32e2ac2d 100644 --- a/src/SimpleSDMLayers.jl +++ b/src/SimpleSDMLayers.jl @@ -23,12 +23,14 @@ include(joinpath("datasets", "sources.jl")) include(joinpath("datasets", "download_layer.jl")) export EarthEnv, WorldClim, BioClim -include(joinpath("datasets", "geotiff.jl")) include(joinpath("datasets", "raster.jl")) +include(joinpath("datasets", "ascii.jl")) +include(joinpath("datasets", "geotiff.jl")) +export geotiff + include(joinpath("datasets", "worldclim.jl")) include(joinpath("datasets", "chelsa.jl")) include(joinpath("datasets", "landcover.jl")) -include(joinpath("datasets", "ascii.jl")) export worldclim, bioclim, landcover include(joinpath("operations", "coarsen.jl")) diff --git a/src/datasets/ascii.jl b/src/datasets/ascii.jl index 917e59be..c792f229 100644 --- a/src/datasets/ascii.jl +++ b/src/datasets/ascii.jl @@ -32,7 +32,7 @@ function ascii(file::AbstractString, datatype::Type{T}=Float64) where {T <: Numb grid[i] = nothing end end - return SimpleSDMPredictor(grid, xl, xl+cs*2ncols, yl, yl+cs*2nrows) + return SimpleSDMPredictor(grid, xl, xl+cs*ncols, yl, yl+cs*nrows) end """ @@ -47,7 +47,7 @@ function ascii(layer::SimpleSDMPredictor{T}, file::AbstractString; nodata::T=con open(file, "w") do io write(io, "ncols $(size(layer, 2))\n") write(io, "nrows $(size(layer, 1))\n") - write(io, "cellsize $(stride(layer)[1])\n") + write(io, "cellsize $(2*stride(layer)[1])\n") write(io, "xllcorner $(layer.left)\n") write(io, "yllcorner $(layer.bottom)\n") write(io, "nodata_value $(nodata)\n") diff --git a/src/datasets/geotiff.jl b/src/datasets/geotiff.jl index 8af0056d..058e3156 100644 --- a/src/datasets/geotiff.jl +++ b/src/datasets/geotiff.jl @@ -76,3 +76,47 @@ function geotiff( return LT(buffer, left_pos-0.5lon_stride, right_pos+0.5lon_stride, bottom_pos-0.5lat_stride, top_pos+0.5lat_stride) end + +function geotiff(layer::SimpleSDMPredictor{T}, file::AbstractString; nodata::T=convert(T, -9999)) where {T <: Number} + array = layer.grid + replace!(array, nothing => NaN) + array = convert(Matrix{T}, array) + dtype = eltype(array) + array_t = reverse(permutedims(array, [2, 1]); dims=2) + width, height = size(array_t) + + # Geotransform + gt = zeros(Float64, 6) + gt[1] = layer.left + gt[2] = 2stride(layer, 1) + gt[3] = 0.0 + gt[4] = layer.top + gt[5] = 0.0 + gt[6] = -2stride(layer, 2) + + # Write + prefix = first(split(last(splitpath(file)), '.')) + ArchGDAL.create(prefix, + driver=ArchGDAL.getdriver("MEM"), + width=width, height=height, + nbands=1, dtype=T, + options=["COMPRESS=LZW"]) do dataset + + band = ArchGDAL.getband(dataset, 1) + + # Write data to band + ArchGDAL.write!(band, array_t) + + # Write nodata and projection info + ArchGDAL.setnodatavalue!(band, nodata) + ArchGDAL.setgeotransform!(dataset, gt) + ArchGDAL.setproj!(dataset, "EPSG:4326") + + # Write ! + ArchGDAL.write(dataset, file, driver=ArchGDAL.getdriver("GTiff"), options=["COMPRESS=LZW"]) + end +end + +function geotiff(layer::SimpleSDMResponse{T}, file::AbstractString; nodata::T=convert(T, -9999)) where {T <: Number} + geotiff(convert(SimpleSDMPredictor, layer), file; nodata=nodata) +end \ No newline at end of file diff --git a/test/ascii.jl b/test/ascii.jl index 2442bbff..7787dffa 100644 --- a/test/ascii.jl +++ b/test/ascii.jl @@ -20,6 +20,10 @@ U = SimpleSDMLayers.ascii("test.asc") @test S.top == U.top @test size(S) == size(U) +geotiff(U, "test.tiff") +@test isfile("test.tiff") + rm("test.asc") +rm("test.tiff") end \ No newline at end of file