Hi Thiago,
Something like this should work:
s <- stack(lapply(1:5, function(x) {
r <- raster(nrows=22, ncols=20, xmn=-58, xmx=-48, ymn=-33, ymx=-22)
r[] <- round(runif(22 * 20, min=0, max=600), digits=0)
r
}))
intervals <- seq(0, 600, 100)
vals <- sapply(1:nlayers(s), function(x) { # x
Hi Thiago,
Done in haste, but I think this might do it (it’s on an 8X8 problem though):
stlist <- lapply(1:8, function(x) {
rl <- stack(lapply(1:8, function(y) {
r <- raster(nrow = 10, ncol = 10)
r[] <- sample(1:100, size = ncell(r), replace = TRUE)
r
}))
rl
})
Hi Milan,
Two possible workarounds I have found that can sometimes clear up a
problem with unclean polygons are to first run either:
gBuffer(your_shape, width = 0)
gUnaryUnion (if you don't mind losing the identities of individual polygons)
And then try run the intersection on the output. It
I am not sure about the mismatch issue, but I thinking merging the
data slot of spatialPolygonsDataFrame with a data frame produces
undesirable results.
I wrote a function a while back that does the merge in such a way that
the problems are avoided, and perhaps this might help. I think there
are
Hi Sam,
One thing you might check is the column names in your csv file. If you
can change them to something simpler, without periods in them, etc,
that might get rid of the problem. There are some threads on this
problem as well, I believe.
Cheers, Lyndon
On Thu, Jul 31, 2014 at 11:27 AM,
Hi Katrina,
I might have misunderstood what you are trying to do, but its seems
the main problem is that you are using the calc function rather than
overlay. Also, I think you don't need the separate function for what
you are trying to achieve, if you are just trying to apply this to two
rasters:
Dear List,
I am trying to create some maps with postscript, which start with a
base map (of Africa, with grey fill, or with a raster with a single
grey color for all non-NA pixels), and then have a raster that
partially covers the area plotted using a divergent color scheme. The
grey background
Hi Ervan,
You can try first do a gbuffer with a width of 0 (can't test now, sorry)
and see if that helps, but for a more effective solution you might want to
try build and call to pprepair, which you can find on github (Barry
Rowlinson pointed me to this very helpful software).
Hope this helps.
Hi Megan,
You might be able to do both the resampling and reprojecting quite
handily with gdal, using a system call out of R.
Here are some lines from one of my scripts (I haven't looked at this
in a while and haven't verified now, so caveat emptor) that might
help:
# lc is a landsat landcover
Hi Megan,
Can you give some steps on how you prepared the data before stacking?
Was there any resampling, reprojecting, etc. involved in your forest
cover layer? I have had similar results with rasters in the past
where NA values were represented by negative numbers. What might have
happened is
Hi Barry,
Thanks very much for the tip. I will try and see if I can implement
that using system calls from my R script, and will report back.
Thanks, Lyndon
On Fri, Dec 14, 2012 at 7:04 PM, Barry Rowlingson
b.rowling...@lancaster.ac.uk wrote:
On Fri, Dec 14, 2012 at 8:16 PM, Lyndon Estes les
Dear List,
I am working on a project to have users maps agricultural fields
(mappingafrica.princeton.edu, if anyone is interested or has comments
on what we are trying to do (will be highly appreciated)), and it
makes uses of rgeos functions quite heavily to perform some quality
control checks on
/ passwd here:
http://cran.r-project.org/web/packages/spacetime/vignettes/stpg.pdf (in
this example, readOGR was not used because only the attribute tables
resulting from spatial and/or temporal selections in postGIS were retrieved)
On 10/10/2012 02:01 PM, Lyndon Estes wrote:
Hi to everyone,
Thanks
Hi to everyone,
Thanks to all of you for the advice on this, which has solved my problem
quite nicely now. Based on your suggestions, I figured out the following:
ogr2ogr was broken for me because I had not completely removed a prior
custom install of gdal that I had set up for hdf4 support,
Dear List,
I was hoping someone could give me some pointers please for trying to read
postGIS polygon geometries into R spatialPolygons.
I unfortunately am not sure how to give a toy example with a postGIS
database, but what I want to do is:
From R, I want to query a postGIS database in order
Hi Robert,
Just wanted to report back on the solution you gave me, which worked
perfectly as you wrote it.
library(raster)
library(rgdal)
r - raster(ncol=10, nrow=10) # EP MODIS pixels
eplc - raster()
set.seed(0)
eplc[] - round( runif(ncell(eplc)) + 3)
eplc[1:11000] - NA
bs -
for how to write a speedier code for addressing this
particular application.
Thanks, Lyndon
error is due to the fact that I ended up with a line for
On Sat, Aug 18, 2012 at 6:38 PM, Lyndon Estes lyndon.es...@gmail.comwrote:
Dear List,
I am currently trying to figure out the the proportions
- writeStop(rcl3)
rcl4 - writeStop(rcl4)
rtot - writeStop(rtot)
plot(stack(rcl3, rcl4, rtot))
On Fri, Aug 24, 2012 at 12:47 PM, Lyndon Estes lyndon.es...@gmail.comwrote:
Hi,
To follow up on my last with an update. I modified the code a bit to run
with my data and set it on its way on our
Dear List,
I am currently trying to figure out the the proportions of two
different landcover classes (classified from Landsat imagery) that
fall within the extent of coarser MODIS pixels (I am trying to isolate
MODIS pixels associated with fairly pure cover types). I want to do
this without
To add to that, the following link has some methods for that. I was doing
the same thing you were doing when I posted for help.
http://r-sig-geo.2731867.n2.nabble.com/Efficient-way-to-obtain-gridded-count-of-overlapping-polygons-td6034590.html
Cheers, Lyndon
On Mon, Jul 2, 2012 at 8:25 AM,
Hi Gavin,
I have been working with some code that does something very similar. In my
case, I am intersecting points with grids of varying resolutions in order
to collect the grid identifier (and coordinates for grid cell centers),
then averaging the results of several variables held by the
Dear List,
I am trying to plot an RGB image using plotRGB, and I was wondering if
there is anyway to turn the NA regions of the image to a different color
(say white or grey)?
For instance, in the following code, how to a turn the first and last 5
columns of the image something other than black?
:13 AM, Lyndon Estes lyndon.es...@gmail.comwrote:
Dear List,
I am trying to plot an RGB image using plotRGB, and I was wondering if
there is anyway to turn the NA regions of the image to a different color
(say white or grey)?
For instance, in the following code, how to a turn the first
It looks like you have not put set up lapply incorrectly.
hadcm3_40-lapply(list.files(
path=
Original_data/AR4_climate_change_ASCII_grids/,pattern=^tas_2040.*a1b.*hadcm3),raster,
crs=prj_latlong)
Try this instead:
fnms - list.files(path= Original_data/AR4_climate_change_ASCII_grids/,
:
Hi
Attaching 2 files.
Thanks
** **
Kanika
** **
*From:* Lyndon Estes [mailto:lyndon.es...@gmail.com]
*Sent:* Friday, 2 March 2012 4:11 p.m.
*To:* Kanika Jhunjhnuwala
*Cc:* r-sig-geo@r-project.org
*Subject:* Re: [R-sig-Geo] reading a list of text files as rasters
Something like this should do the trick
plot(shp) # Plot your shapefile in Lambert projection
plot(spTransform(pts, shp@proj4string), add = T) # Add points converted to
Lambert
Cheers, Lyndon
On Sun, Feb 26, 2012 at 12:24 PM, gregor.hochsch...@gmx.de wrote:
Hi,
I am not sure whether I
Hi Johannes,
For number two, you could also use this:
?dropLayer
Cheers, Lyndon
On Tue, Feb 14, 2012 at 7:51 AM, Johannes Radinger jradin...@gmx.at wrote:
Dear R-sig-geo List!
I have two little questions concerning RasterBricks from the Raster
package.
1) I have a list of several
Hi Adam,
Here is code for something very similar that I was doing, which should help
with what you want to do. In this case, I made grids of varying resolution
to cover the extent of the area defined by my points, and then collected
the grid ids for each resolution and bound them to the
Hi Ana,
Based on a previous post I found, it looks like the package SDMTools could
do what you want via the PatchStat function:
http://cran.r-project.org/web/packages/SDMTools/index.html
As to the question put regarding why use R for GIS work, I can give my
perspective as someone who has
meters) the polygon should correspond to that polygon. So, it this
way I could change the precision of the estimation, don't you?
Than you in advance!
José
--
*De:* Lyndon Estes lyndon.es...@gmail.com
*Para:* Jose Bustos Melo jbustosm...@yahoo.es
*CC:* r-sig-geo
- stackSelect(rb2, ind.r)
I.e., for each cell, use the value in ind.r to select the layer in rb2 from
which to take the value for the output layer.
Best, Robert
On Thu, Dec 1, 2011 at 5:55 PM, Lyndon Estes lyndon.es...@gmail.com wrote:
Dear List,
I was hoping to get some advice on the best way
21, 2011 at 8:40 AM, Lyndon Estes lyndon.es...@gmail.com
wrote:
Hi Swen,
Thanks for the suggestion. That brings to mind an earlier problem I
had with a different functions substituting values on large rasters,
and I got around it by writing the following lines into the code:
if(raster
Hi Swen,
Thanks for the suggestion. That brings to mind an earlier problem I
had with a different functions substituting values on large rasters,
and I got around it by writing the following lines into the code:
if(raster:::.toDisk() != TRUE) {
setOptions(todisk = TRUE)
cat(We don't
Dear List,
I am having trouble creating a brick out of a raster stack with 9
grids of 2.5 mb each. My goal in creating the brick is to hopefully
speed up calculation of quantiles for each grid location, which runs
very slowly as follows:
q.10 - calc(stack1, function(x) quantile(x, 0.1))
I was
Hi Maria,
Please see my inserts below.
---
The script is as following:
install.packages(sp)
library(sp)
#load the map #
mapchina- url(http://www.gadm.org/data/rda/CHN_adm1.RData;)
load(mapchina)
mapchina=gadm
gadm$NAME_1 #shows the name of the Chinese regions
#create the vector
Hi Sam,
I wrote the following function a while ago to account for this problem
of dataframe row reordering relative to the polygon IDs:
joinAttributeTable - function(x, y, xcol, ycol) {
# Merges data frame to SpatialPolygonsDataFrame, keeping the correct
order. Code from suggestions at:
#
, at 12:27 PM, Lyndon Estes wrote:
Hi Sam,
I wrote the following function a while ago to account for this problem
of dataframe row reordering relative to the polygon IDs:
joinAttributeTable - function(x, y, xcol, ycol) {
# Merges data frame to SpatialPolygonsDataFrame, keeping the correct
order
Hi Jonathan,
You can do this as follows:
x - yourSpatialPolygons
xl - as(x, SpatialLines)
e.g.
library(maptools)
library(rgdal)
nc - readShapePoly(system.file(shapes/sids.shp,
package=maptools)[1], proj4string=CRS(+proj=longlat +datum=NAD27))
quartz()
plot(nc[1, ], col = blue)
nc1l - as(nc[1,
Hi Sara,
The following is slightly different from what you are trying to do,
but should give you an area value for cells of a given value, without
having to stack your rasters:
r - raster(ncol = 10, nrow = 10)
res(r) - 10 # Giving it a square grid cell, which remakes raster
into 18 rows 36
Hi Robert,
Many thanks for the answer. I tried this again using median rather
than Median, and got the following results on my example:
setOptions(todisk = TRUE) # I used this option to rule out memory issues
system.time(med - calc(stk, median, filename = med.tif, datatype =
INT2S, overwrite =
Dear List,
I am trying to extract the per-pixel median from a raster stack, and
am finding that the process takes a very long time relative to taking
the per-pixel mean. It is long enough that I have never actually seen
whether it successfully completes on my data, but this dummy example
shows
: mdsum...@gmail.com
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
Lyndon Estes
Research Associate
Woodrow Wilson School
Princeton University
+1-609-258-2392 (o)
+1-609-258-6082 (f)
+1-202-431
Hi Colin,
How many colors are actually displayed? I find that the legend/scale
bar often shows a range of values, even though there is only one value
in the raster.
For instance if I do this:
dummy - raster(extent(-100, 100, -200, 200))
dummy[] - 1
plot(dummy)
I have the same values as you in
others bring data.
[[alternative HTML version deleted]]
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
Lyndon Estes
Research Associate
Woodrow Wilson School
Princeton University
+1
Hi Robert,
My apologies for an even longer delay this time around. I have gotten
back to this now and followed your suggestions. Responses are
interspersed in (shortened) text of previous messages.
That is *horrible*. I am not sure what is going on here. This is my
rather empirical test to
/listinfo/r-sig-geo
--
Lyndon Estes
Research Associate
Woodrow Wilson School
Princeton University
+1-609-258-2392 (o)
+1-609-258-6082 (f)
+1-202-431-0496 (m)
les...@princeton.edu
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch
://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
Lyndon Estes
Research Associate
Woodrow Wilson School
Princeton University
+1-609-258-2392 (o)
+1-609-258-6082 (f)
+1-202-431-0496 (m)
les...@princeton.edu
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https
to the area covered by
sum.tab.c works fine), I don't want to burden anyone with such a large
file.
Thanks in advance for any advice you can provide.
Cheers, Lyndon
--
Lyndon Estes
Research Associate
Woodrow Wilson School
Princeton University
+1-609-258-2392 (o)
+1-609-258-6082 (f)
+1-202-431-0496 (m
Hello,
I am interested in creating a raster that summarizes the number of
overlapping polygons that underlie each cell (in order to display how
many species occur at that point).
I have come up with a function that seems to work, but I was wondering
if I had missed one that already exists, or if
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
Lyndon Estes
Research Associate
Woodrow Wilson School
Princeton University
+1-609-258-2392 (o)
+1-609-258-6082 (f)
+1-202-431-0496 (m)
les...@princeton.edu
___
R-sig-Geo
= -35, h = 10, w = 10, size = 2.5)
GridTopology and related help pages don't mention how to use degrees
instead of metres (or my R immaturity prevented me from finding)...
Any tips to help me working on that adaptation?
Best regards,
Thiago.
--- On Tue, 15/2/11, Lyndon Estes les
51 matches
Mail list logo