Geoprocessing functionalities in Julia


#1

One more general question - in case some of the other users come across some some helpful packages and webpages giving some advice for geoprocessing code-examples in JULIA, it would be great if you could share them here.

I am currently getting acquainted with the ArchGDAL package which provides most of the basic GDAL functionalities in a quite comfortable wrapper.

https://github.com/yeesian/ArchGDAL.jl

However, I am still searching for many base-funtions like for example extracting raster values located within a Polygon etc. (in R this would be the extract()-function of the raster package)

I will for now probably simply switch to R for these steps, but in case there are solutions for JULIA available, it would of course be nice to implement everything in one environment.


#2

I think you will not find GDAL wrappers which are as convenient as R raster yet in Julia. However, to solve your problem, you could try to use this rasterize function I once made in Julia using GDAL:

# GDALRasterize
function rasterize!(outar,shapefile;x0=-180.0,y0=-90.0,x1=180.0,y1=90.0)
    sx,sy = size(outar)
    ds_point = GDAL.openex(shapefile, GDAL.GDAL_OF_VECTOR, C_NULL, C_NULL, C_NULL)
    ds_point == C_NULL && error("Could not open shapefile")
    options = GDAL.rasterizeoptionsnew(["-of","MEM","-ts",string(sx),string(sy),"-te",string(x0),string(y0),string(x1),string(y1),"-a_srs","WGS84"], C_NULL)
    ds_rasterize = GDAL.rasterize("data/point-rasterize.mem", Ptr{GDAL.GDALDatasetH}(C_NULL), ds_point, options, C_NULL)
    GDAL.rasterizeoptionsfree(options)
    GDAL.close(ds_point)

    band = GDAL.getrasterband(ds_rasterize,1)

    xsize = GDAL.getrasterbandxsize(band)
    ysize = GDAL.getrasterbandysize(band)
    dtype = GDAL.getdatatypename(GDAL.getrasterdatatype(band))
    @assert (xsize,ysize)==(sx,sy)
    GDAL.rasterio(band, GDAL.GF_Read, 0, 0, xsize, ysize,
    outar, xsize, ysize, GDAL.GDT_Float64, 0, 0)
    GDAL.close(ds_rasterize)
end

This will rasterize a given shapefile (path) into the given array outar, where the bounding box of the output array can be specified. You can use this resulting array as a mask to extract the values inside the shapefile from the ESDC. Let me know if you need further help with this approach, I could certainly make up a more detailed example.