You can use raster package to download WorldClim data, see ?getdata to know about resolution, variables and coordinates.
r <- getData("worldclim",var="bio",res=10)
Bio 1 and Bio12 are mean anual temperature and anual precipitation:
r <- r[[c(1,12)]]
names(r) <- c("Temp","Prec")
I create random points as ...
Use calc() to apply functions over a raster object, such as Raster, RasterStack or RasterBrick:
mean <- calc(STACK1, fun = mean)
If you have na values in cells, add na.rm =T:
mean <- calc(STACK1, fun = mean, na.rm = T)
Alternatively, you can use stackApply also:
mean <- stackApply(STACK1, indices = rep(1,nlayers(STACK1)), fun = "mean", na.rm = ...
You need to use the beginCluster and endCluster functions of the raster package. See the example below.
# Make test data
r <- raster(ncol=36, nrow=18)
r <- 1:ncell(r)
s <- stack(r, sqrt(r), r/r)
cds1 <- rbind(c(-180,-20), c(-160,5), c(-60, 0), c(-160,-60), c(-180,-20))
cds2 <- ...
Convert df to a sp object and use sp = TRUE as argument.
What you have:
r <- raster(ncol=36, nrow=18)
r <- 1:ncell(r)
xy <- cbind(-50, seq(-80, 80, by=20))
##  626 554 482 410 338 266 194 122 50
What you need:
extract(r, SpatialPoints(xy), sp = T)
## class : SpatialPointsDataFrame
## features : 9
This is not related to R 4.0.1 but to rgdal 1.5-8 and the migration to gdal 3 and proj 6. This is a very long and complex process that impact hundreds maybe thouthands of packages. All the packages are not yet up-to-date with what is coming.
You can have a look to ?rgdal::set_thin_PROJ6_warnings() to eliminate those warnings.
Edit after Roger's comment: See ...
The raster::extract function, when applied to polygons or with the buffer argument, returns a list object where each element in the list contains a vector of the raster values intersecting the polygon. If the input raster object was a stack or brick, containing multiple rasters, the list elements are a matrix rather than a vector.
Providing the fun argument ...
I get much faster results with velox if I crop the raster before running extract, e.g.:
r <- velox("testras_so.tif")
I've also been working on a package with an optimized extract function that may be of interest:
exact_extract(ras, poly) # get a matrix with weights and values,
You can do Layer... Save As... and choose the CSV output format. Choose Geometry type as POINT and Geometry as AS_XY:
Then you'll get a CSV like:
You need to transform your points to the coordinate system of the raster. You could warp the raster to your points' coordinate system but warping rasters is a bit messy.
First create a SpatialPoints object from your coordinates and tell R it is in lat-long coordinates (4326):
pts2 = SpatialPoints(points,...
3D analyst's Add Surface Information will add a Z field to your vector data with the data value from the overlapping raster layer:
Interpolates surface elevation properties for point, multipoint, and
That's for v10, I didn't catch which version of ArcGIS you were using.
From what I can find there doesn't appear to be an existing solution for this exact situation, but I still wanted to be able to do this in QGIS, so I took the plunge into python scripting.
A guide for writing processing algorithms can be found here https://docs.qgis.org/2.18/en/docs/user_manual/processing/scripts.html
To use this code open up the ...
Take a look at ee.FeatureCollection.aggregate_array:
Aggregates over a given property of the objects in a collection,
calculating a list of all the values of the selected property.
should bring you an array of "means".
To provide the response based on @artwork21 comment.
To identify areas in one shapefile that contain the areas in a different shapefile as related to your soils data. You can than extract the data in your attribute table using a Spatial Join feature. Here is a tutorial using QGIS that you can reference Performing Spatial Joins
I've included a few of the ...
+proj=longlat is in angular units, degrees by default.
Better if you can project your dataset before, to a conformal projection (to preserve the circular buffer shape) appropriate for your work area (so that the distance deforms as little as possible).
Run "Select by location" from your processing toolbox and choose your polygonlayer as "Select features from", your linelayer as "By comparing to the features from" and intersect as "geometric predicate". Your layers do not need to be in the same CRS, that is not an issue.
If you want to extract these features. You can ...
You can create a rectangle with the boundary box of your polygon. I use rgdal package to make a reproducible example, but you could use only raster package with shapefile() function instead of readOGR():
dsn <- system.file("vectors", package = "rgdal")
polygon <- readOGR(dsn=dsn, layer="scot_BNG")
The coordinate system information is incorrect. You can usually find the information in a .e00 file. Here's what it says:
spheroid = GRS80
central meridian/longitude of origin = -96.0
standard parallel 1 = 29.5
standard parallel 2 = 45.5
latitude of origin = 23.0
false easting/false northing = 0.0
so the PROJ.4 string should be:
+proj=aea +lat_1=29.5 +...
You can compress it further by right-clicking the clipped tiff in the layer list and selecting Save As...
In that window, the format at the top should be kept to GTiff and set the CRS to whatever one you are using.
Go down to Create Options and change the Profile to High Compression.
If you want to do this in one step, try using the Clipper tool for ...
I ended up using an approach based on the snowfall package. It is quite simple, works really good and the point extraction function is as fast as the number of cores that you can use. The approach I used was inspired by this post, and here is my reproducible example:
# Create date sequence
idx <- seq(as.Date("2010/1/1")...
If you are interested in the z-coordinate of each point, you should convert the line into points. You can use for example the "Extract Nodes" tool in QGIS for this task. After that open the Field Calculator for the point layer, paste z($geometry) like in the image below and write the result in a field in the attribute table.
I think this should work:
var dates = ee.List(["property1", "property2"]);
f = ee.Feature(f);
return ee.Feature(null, f.toDictionary(dates));
First of all, you have some syntax errors in your example code. You need to give the inputs of the function, here that is "(f)", after the "function" statement. Then ...
Most of the tool have a description on the right side of the window (if not visible click on the arrow in the upper right corner)
In your case the description read :
This algorithm takes a line or polygon layer and generates a point layer with points representing the vertices in the
input lines or polygons. The attributes associated to each ...
There is a possibility of using a 'Virtual Layer' through Layer > Add Layer > Add/Edit Virtual Layer...
Let's assume there two layers, a line layer 'lines' and a polygon layer 'polygons'.
With the following query, it is possible to extract polygons from a layer crossed by a lines layer.
FROM "polygons" AS p, "lines" AS l
You can do this with Hawthorne Beyer's free Geospatial Modelling Environment, (GME, formerly known as Hawth's Tools). There is a tool in there, Intersect Points With Raster, which as its name implies, acts like the Intersect tool in ArcGIS but allows you to intersect a point layer with a raster, like the Extract Values to Points tool. You can also apply an ...
I have been facing the same issue some time ago. For me, its not the 'Select by attribute' algorithm, but the 'Extract by attribute', that delivers the desired result. I think the reason is, that 'Select by attribute' only creates a selection in QGIS, but does not return any data. 'Extract by attribute' creates a new dataset, which can be used for further ...