joshobrien / gdalutilities Goto Github PK
View Code? Open in Web Editor NEWHome Page: https://joshobrien.github.io/gdalUtilities/
Home Page: https://joshobrien.github.io/gdalUtilities/
This package depends on (depends, imports or suggests) raster and one or more of the retiring packages rgdal, rgeos or maptools (https://r-spatial.org/r/2022/04/12/evolution.html, https://r-spatial.org/r/2022/12/14/evolution2.html). Since raster 3.6.3
, all use of external FOSS library functionality has been transferred to terra, making the retiring packages very likely redundant. It would help greatly if you could remove dependencies on the retiring packages as soon as possible.
I want to wrap the gdal functions in a function with some special default parameters, like below:
myfun<- function(input,out="out.tif",projwin=NULL,...){
if(is.null(projwin) ) projwin=substitute()
gdalUtilities::gdal_translate(input, dst_dataset=out, projwin=projwin,...)
}
but it does not work with an error : Error in FUN(X[[i]], ...) : argument "val" is missing, with no default
I can work aroud by:
if(missing(projwin)){
gdalUtilities::gdal_translate(fname,dst_dataset=out,...)
}else{
gdalUtilities::gdal_translate(fname,dst_dataset=out,projwin=projwin)
}
But it would be not efficient if there is more parameters for cycling.
Would you please tell me what is the right way?
When I run gdalwarp with function = "sum" it creates a single column at the easternedge of the output raster with very large values. It seems that the values in that column are equal the sum of all the pixels at that latitude. I do not get this error when I run gdalwarp from the OSGeo4W command line. It's possible gdalUtilities is using an older version of gdal, but I don't know how to check this.
Here's my command:
gdalUtilities::gdalwarp(srcfile = in.tif", dstfile = "out.tif", s_srs = "EPSG:4326", t_srs = "EPSG:4326", tr= c(0.1,0.1), te = c(-180,-56,180,84), r = "sum", overwrite = TRUE)
And here is the result:
I am using a gridded population product: https://sedac.ciesin.columbia.edu/data/set/popdynamics-1-km-downscaled-pop-base-year-projection-ssp-2000-2100-rev01
gdalUtilities::gdal_rasterize(shp, f <- tempfile(fileext='.tif'),
tr=c(0.01, 0.01), a='ID_2')
## OGR data source with driver: ESRI Shapefile
## Source: "/Library/Frameworks/R.framework/Versions/3.6/Resources/library/raster/external/lux.shp", layer: "lux"
## with 12 features
## It has 5 fields
## Error in CPL_gdalrasterize(source, destination, options) :
## cannot write to destination dataset
## In addition: Warning message:
## In CPL_gdalrasterize(source, destination, options) :
## GDAL Error 4: /var/folders/lj/k338sqgj0c5cxzftty8mhbn0f68cd3/T//RtmpcUEXeB/file39bb11f86012.tif: No such file or directory
Equivalent with gdalUtils
does create the output tif:
gdalUtils::gdal_rasterize(shp, f <- tempfile(fileext='.tif'),
tr=c(0.01, 0.01), a='ID_2')
I believe output raster creation should be possible since GDAL 1.8.0.
Hello, we depend on your package and wanted to check if you plan to release an updated version for CRAN soon.
Passing multiple options to co
(for functions that accept it) seems to be ineffective, with only the first option respected.
e.g., compare:
gdalUtilities::gdal_translate(
system.file(package='gdalUtilities', 'extdata', 'tahoe.tif'),
f <- tempfile(fileext='.tif'), co=c('TILED=YES', 'COMPRESS=LZW'))
gdalUtilities::gdalinfo(f)
## .......
## Image Structure Metadata:
## INTERLEAVE=PIXEL
## Corner Coordinates:
## Upper Left (-119.9328072, 39.2914138) (119d55'58.11"W, 39d17'29.09"N)
## Lower Left (-119.9328072, 39.2892247) (119d55'58.11"W, 39d17'21.21"N)
## Upper Right (-119.9306180, 39.2914138) (119d55'50.22"W, 39d17'29.09"N)
## Lower Right (-119.9306180, 39.2892247) (119d55'50.22"W, 39d17'21.21"N)
## Center (-119.9317126, 39.2903192) (119d55'54.17"W, 39d17'25.15"N)
## Band 1 Block=256x256 Type=Byte, ColorInterp=Red
## Band 2 Block=256x256 Type=Byte, ColorInterp=Green
## Band 3 Block=256x256 Type=Byte, ColorInterp=Blue
and
gdalUtilities::gdal_translate(
system.file(package='gdalUtilities', 'extdata', 'tahoe.tif'),
f2 <- tempfile(fileext='.tif'), co=c('COMPRESS=LZW', 'TILED=YES'))
gdalUtilities::gdalinfo(f2)
## .......
## Image Structure Metadata:
## COMPRESSION=LZW
## INTERLEAVE=PIXEL
## Corner Coordinates:
## Upper Left (-119.9328072, 39.2914138) (119d55'58.11"W, 39d17'29.09"N)
## Lower Left (-119.9328072, 39.2892247) (119d55'58.11"W, 39d17'21.21"N)
## Upper Right (-119.9306180, 39.2914138) (119d55'50.22"W, 39d17'29.09"N)
## Lower Right (-119.9306180, 39.2892247) (119d55'50.22"W, 39d17'21.21"N)
## Center (-119.9317126, 39.2903192) (119d55'54.17"W, 39d17'25.15"N)
## Band 1 Block=400x6 Type=Byte, ColorInterp=Red
## Band 2 Block=400x6 Type=Byte, ColorInterp=Green
## Band 3 Block=400x6 Type=Byte, ColorInterp=Blue
dryrun
shows that the required (afaik) -co
are not repeated for each creation option.
e.g.
gdalUtilities::gdal_translate(
system.file(package='gdalUtilities', 'extdata', 'tahoe.tif'),
tempfile(fileext='.tif'), co=c('COMPRESS=LZW', 'TILED=YES'),
dryrun=TRUE)
## [1] "gdal_translate /Library/Frameworks/R.framework/Versions/3.6/Resources/library/gdalUtilities/extdata/tahoe.tif /var/folders/lj/k338sqgj0c5cxzftty8mhbn0f68cd3/T//RtmpJjiLbM/fileb30233f73487.tif -co COMPRESS=LZW TILED=YES"
Tests are not found because the necessary testthat.R
file is forgotten:
---> Building R-gdalUtilities
---> Testing R-gdalUtilities
* using log directory ‘/opt/local/var/macports/build/_opt_PPCSnowLeopardPorts_R_R-gdalUtilities/R-gdalUtilities/work/gdalUtilities/gdalUtilities.Rcheck’
* using R version 4.3.1 (2023-06-16)
* using platform: powerpc-apple-darwin10.0.0d2 (32-bit)
* R was compiled by
gcc-mp-12 (MacPorts gcc12 12.3.0_0) 12.3.0
GNU Fortran (MacPorts gcc12 12.3.0_0) 12.3.0
* running under: OS X Snow Leopard 10.6
* using session charset: UTF-8
* checking for file ‘gdalUtilities/DESCRIPTION’ ... OK
* checking extension type ... Package
* this is package ‘gdalUtilities’ version ‘1.2.5’
* package encoding: UTF-8
* checking package namespace information ... OK
* checking package dependencies ... OK
* checking if this is a source package ... OK
* checking if there is a namespace ... OK
* checking for executable files ... OK
* checking for hidden files and directories ... OK
* checking for portable file names ... OK
* checking for sufficient/correct file permissions ... OK
* checking whether package ‘gdalUtilities’ can be installed ... OK
* checking installed package size ... OK
* checking package directory ... OK
* checking DESCRIPTION meta-information ... OK
* checking top-level files ... OK
* checking for left-over files ... OK
* checking index information ... OK
* checking package subdirectories ... OK
* checking R files for non-ASCII characters ... OK
* checking R files for syntax errors ... OK
* checking whether the package can be loaded ... OK
* checking whether the package can be loaded with stated dependencies ... OK
* checking whether the package can be unloaded cleanly ... OK
* checking whether the namespace can be loaded with stated dependencies ... OK
* checking whether the namespace can be unloaded cleanly ... OK
* checking loading without being on the library search path ... OK
* checking dependencies in R code ... OK
* checking S3 generic/method consistency ... OK
* checking replacement functions ... OK
* checking foreign function calls ... OK
* checking R code for possible problems ... OK
* checking Rd files ... OK
* checking Rd metadata ... OK
* checking Rd cross-references ... OK
* checking for missing documentation entries ... OK
* checking for code/documentation mismatches ... OK
* checking Rd \usage sections ... OK
* checking Rd contents ... OK
* checking for unstated dependencies in examples ... OK
* checking examples ... OK
* checking PDF version of manual ... OK
* DONE
Status: OK
Reproducible example:
m <- matrix(rnorm(n=10^4,m=0,sd=2),nrow=100)
m[1:10,] <- NA
imatest <- raster(m)
cellStats(imatest,range)
imatest <- brick(imatest,imatest*2,imatest+10)
extent(imatest) <- c(0,100,0,100)
projection(imatest) <- CRS("+init=epsg:25831")
imatest
class : RasterBrick
dimensions : 100, 100, 10000, 3 (nrow, ncol, ncell, nlayers)
resolution : 1, 1 (x, y)
extent : 0, 100, 0, 100 (xmin, xmax, ymin, ymax)
coord. ref. : +init=epsg:25831 +proj=utm +zone=31 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
data source : in memory
names : layer.1, layer.2, layer.3
min values : -7.844195, -15.688390, 2.155805
max values : 7.681265, 15.362530, 17.681265
writeRaster(imatest,filename="imatest",format="GTiff",NAflag=-9999,overwrite=TRUE)
qt <- quantile(imatest,probs=c(0.02,0.98),na.rm=TRUE)
qt
scal <- 1
for(i in 1:nlayers(imatest)){
scal <- c(scal,c(qt[i,],c(1,255)))
}
scal <- scal[-1]
> scal
2% 98% 2% 98%
-4.133383 4.107563 1.000000 255.000000 -8.266766 8.215126 1.000000 255.000000
2% 98%
5.866617 14.107563 1.000000 255.000000
> gdal_translate(src_dataset="imatest.tif",
+ dst_dataset="imateststretchgdal2.tif",
+ ot="Byte",scale=scal,dryrun=TRUE)
Now observe that the scale parameter is wrong:
[1] "gdal_translate imatest.tif imateststretchgdal2.tif -ot Byte -scale -4.13338322257502 4.10756297485811 1 255 -8.26676644515004 8.21512594971622 1 255 5.86661677742498 14.1075629748581 1 255"
scale should be repeating "-scale" for each band:
-scale -4.13338322257502 4.10756297485811 1 255 -scale -8.26676644515004 8.21512594971622 1 255 -scale 5.86661677742498 14.1075629748581 1 255"
and thus the actual run yields an error:
> gdal_translate(src_dataset="imatest.tif",
+ dst_dataset="imateststretchgdal2.tif",
+ ot="Byte",scale=scal)
Warning message:
In CPL_gdaltranslate(source, destination, options) :
GDAL Error 6: Unknown option name '-8.26676644515004'
Same problem if the scale values are entered in matrix form:
> mscal <- matrix(scal,ncol=3, byrow=FALSE)
> mscal
[,1] [,2] [,3]
[1,] -4.133383 -8.266766 5.866617
[2,] 4.107563 8.215126 14.107563
[3,] 1.000000 1.000000 1.000000
[4,] 255.000000 255.000000 255.000000
> gdal_translate(src_dataset="imatest.tif",
+ dst_dataset="imateststretchgdal2.tif",
+ ot="Byte",scale=mscal,dryrun=TRUE)
[1] "gdal_translate imatest.tif imateststretchgdal2.tif -ot Byte -scale -4.13338322257502 4.10756297485811 1 255 -8.26676644515004 8.21512594971622 1 255 5.86661677742498 14.1075629748581 1 255
Thanks for this useful package in prevision of future retirement of rgdal and rgeos.
While using function gdalinfo
, I noticed that it was unable to resolve/expand "~" for user's home directory, at least on Ubuntu 22.04.
Example (tested on Ubuntu 22.04 only):
image <- system.file("extdata", "europe.tif", package="gdalUtilities")
image_copy <- file.path("~", basename(image))
file.copy(image, image_copy)
info <- gdalUtilities::gdalinfo(image_copy)
if(length(info)==0)
warning("It should not be character(0)")
gdalUtilities::gdalinfo(path.expand(image_copy))
# returns the good info
# remove image file copy
file.remove(image_copy)
It would be nice to have the tilde expanded automatically as it is the usual behavior in R.
Also, if file does not exists, gdalinfo just returns a NULL
result. I think it would be expected to send an file not found error.
This issue archives email correspondence from March of 2021 which should have been handled here as an issue. My emails will be marked with a (Josh) and those of the user who contacted me with an (A---).
Hello. I have a multlayer raster and only need first one, or other index, but only one layer.
I tried with several parameters the following function and doesn't work.
gdalUtilities::gdalwarp(srcfile = inraster, dstfile = outraster,
ot = 'Float64',
#dstband = 1, #not working
# b = 1, #not working
# srcband = 1, #not working
tr = rep(max(pixsize), 2),
srcnodata = nd,
dstnodata = -9999)
Thanks
When I try to run gdalwarp
in R (conda environment in VSCode), I get the following. After looking up the error I thought it might be due to available memory or the packages, but the former was fine and for the latter I only had this package plus terra
loaded. I've gotten this error it seems sporadically as well (i.e. a few days ago I had no issue running).
I tried simply using dryrun=TRUE
and then running gdal directly using system()
, which runs with no issue. This makes me think there's something in the package or gdalUtilities::gdalwarp()
function that's interfering, but I'm not sure what it could be.
*** caught segfault ***
address 0xc, cause 'invalid permissions'
Traceback:
1: CPL_gdalwarp(source, destination, options, oo, doo, config_options, quiet, "-overwrite" %in% options)
2: gdal_utils("warp", source = srcfile, destination = dstfile, options = opts, config_options = config_options)
3: gdalwarp(srcfile = f[X], dstfile = outPath, t_srs = "EPSG:32617", te = c(extCommon[1], extCommon[3], extCommon[2], extCommon[4]), tr = c(resMin, resMin), tap = TRUE, r = "cubic", overwrite = TRUE, dryrun = FALSE)
4: FUN(X[[i]], ...)
5: lapply(X = X, FUN = FUN, ...)
6: sapply(1:length(orth$f), standardizeOrthos, f = orth$f, orthos = orth$orthos, pathOriginal, pathStandard, extCommon, resMin)
Possible actions:
1: abort (with core dump, if enabled)
2: normal R exit
3: exit R without saving workspace
4: exit R saving workspace
Selection: 1
R is aborting now ...
zsh: segmentation fault R
A declarative, efficient, and flexible JavaScript library for building user interfaces.
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google ❤️ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.