Tag Info

Hot answers tagged

5

You can: Load required libraries: library(raster) library(rgdal) Read rasters: r1 = raster("./dir/r1.tif") r2 = raster("./dir/r2.tif") Resample to the finer grid r.new = resample(r1, r2, "bilinear") If required (for masking), set extents to match ex = extent(r1) r2 = crop(r2, ex) Removed unrequired data r.new = mask(r.new, r2)


5

What you want is a conditional calculation: return the value of r whenever r and r1 are equal and otherwise set the output to NA. The cell-by-cell arithmetic operations seem to be fastest. (They are much faster than, say, using mask or the reclassification functions.) Since they do not appear to offer an actual conditional operator, use two time-honored ...


3

The official definition for EPSG:3857 is +proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext +no_defs As you see, it is calculated on a sphere (a=b). Your formula is a bit contradictory: You define a lat-long coordinate system on the WGS84 ellispoid, then you add +init=epsg:3857 which should ...


3

There is no real difference other than the class of the returned object. The raster intersect function is a helper function that, for polygons, calls gIntersection from rgeos (not rgdal). I would recommend using raster's intersect functions because it will save you some steps in getting back to a SpatialPolygonsDataFrame object. One good way to explore these ...


2

Looking at the NaturalEarth data page for countries offers you the choice of boundaries with and without lakes which leads me to think you have the one with lakes (but it's hard to tell). For what you are doing you may be better off with the Coastline or Land from the physical set.


2

You can use either rgdal or raster (with an additional step of coersion) to export the prediction, prediction variance or prediction standard deviation as rasters. The sums of squares is stored as a vector where the fit and experimental variogram models are data.frame objects. You will have to attribute a column and output them into a single data.frame or ...


2

I played with assign(), ls(), and mget() to accomplish something that I believe will improve your workflow. First I use ls() to get a list of all environment variables that start with "p": names_poly <- ls(pattern='^p.') I used combn() to find all the unique combinations of polygons combos <- combn(names_poly,2) I looped over combos using get() ...


2

API conditions constantly change, but this should work right now. OSM: devtools::install_github("hrbrmstr/nominatim") library(nominatim) b1 <- osm_geocode("Berlin, Germany") b1[c("lat", "lon")] Yahoo: devtools::install_github("trestletech/rydn") library(rydn) options(RYDN_KEY="yourAPIkey", RYDN_SECRET="yourSecret") b2 <- find_place("Berlin, ...


2

No features at that location! center_bbox takes the arguments center_lon, center_lat, width, height - looks like you're requesting a box centred on 53.3528, -0.4468. That's off the coast of Mogadisho, all blue water, and no features, on OpenStreetMap: Looking at 'bb': left bottom right top 51.3213580 -0.4784529 51.3842420 ...


2

If you reproject a raster with labels, you will obviously get squeezed labels. The only way to avoid this is to render the raster from vector data directly into the desired projection. You might want to look into mapnik, tilemill or maperitive to do this from Openstreetmap raw data (which is vector data). The R openstreetmap package only offers raster ...


1

so, as whuber gave me the nice quick piece of code to create a raster with the values that are the same between rasters, i thought i'd finish the entire job; create rasters and subtract one from the other to get a 'change' raster (full of 0s that need examining as well); r <- raster(ncol=10,nrow=10) r[] <- sample(c(1,2,4,8),size=100,replace=T) ...


1

Have you tried unionSpatialPolygons() from the maptools package?


1

You would calculate it exactly as you would any variogram: by estimating the values of (Z(i) - Z(j))^2/2 as a function of the distance between data locations i and j. You have many choices of distance, but the natural ones would be either distance along the routes or travel time along the routes. If the routes are one-way, additional techniques borrowed ...


1

Convert points to individual lines using Points to Line (Data Management). Place points along the line at regular interval. To do so you might use linear referencing, alternatively search this site, something will pop-up. Calculate distance for each point along the line (Chainage in below table) . Search this site, let me know if fail, I’ll post script or ...


1

Use ceiling. Here is an example in which a raster with values in the interval [1,180] is generated randomly and reclassified. The code then produces a scatterplot of the original and reclassified values to show it has worked correctly. To make it clear that the endpoints of each class are being treated properly (because this is where bugs can arise), the ...


1

This question is similar to: Clip raster by raster with data extraction and resolution change, but coming from a different angle. However, I think the answer is likely the same. First off, choose which raster you wish to be definitive. I'll repeat my previous answer here for ease: Load required libraries: library(raster) library(rgdal) Read rasters: r1 ...


1

Geodetic area calculations are not part of GDAL. You would probably need to write your own function to do it and exact solutions are computationally intensive. You're better off using an Equal-Area projection, such as Albers, Azimuthal, or Lambert azimuthal (https://en.wikipedia.org/wiki/Category:Equal-area_projections). For largish regions in Africa, ...


1

The "landsat-util" Python package used to support querying Landsat images based on a boundary defined by a shapefile until version 0.5.0, where it was removed to get rid of some dependencies. You could try using a version before that, e.g. v0.2.0. It's only a year old, so there is a chance it might still work. The good news is that the package's devs are ...


1

The following code might seem a little long-winded as it represents a manual approach rather than relying on OpenStreetMap, but maybe it's of any help to you anyway. I took the country boundaries and the referring country labels from the wrld_simpl dataset (class 'SpatialPolygonsDataFrame', projected in EPSG:4326) that comes with maptools. The shapefile data ...


1

Just in case anyone has come here for a problem like this one here are a few notes: • epsg: 3577 (as Steven Kay used) is the correct spatial reference in this case (Australia) because it preserves area. Spatial references that don’t exactly preserve area (like epsg:4283) will return a slightly wrong area. Everything will need to be re-projected to the ...


1

from my understanding of your question, and looking at the data, you want to assign a population to each polygon. In pseudocode, something like for each polygon popn <- 0 for each intersection(polygon,pixel outline) fraction <- calculate area of intersection (as percentage of pixel outline value) popn <- popn + ...


1

Thanks everyone. Converting the layer to a raster and exporting as a text file was exactly what I needed for the time being. I used the convert to Raster function in the toolbox. . You can then use the raster to ASCII function to save the data as a txt file which can then be read by my r package.



Only top voted, non community-wiki answers of a minimum length are eligible