[R-sig-Geo] shortestPath in a loop

marta azores martazores at gmail.com
Thu Feb 9 17:38:32 CET 2017


Hi everyone,
the problem was in the projection (+proj=lae) and in the transitionFunction
(mean). If we use the utm projection and the min transitionFunction. It
works properly.
Thanks for your help, in all the process.
:)

Marta

2017-02-03 16:18 GMT-01:00 Chris Reudenbach <reudenbach at uni-marburg.de>:

> Hi Marta,
>
> I am bit of busy wo I will answer later.
>
> cheers chris
>
>
>
> On 03.02.2017 16:32, marta azores wrote:
>
>> Hi Chris,
>> Thanks for your answer. Your solution is perfect. I calculated the length
>> of the line with the old raster (azoTS1), and the results are close to
>> real
>> number.  It's great.(image: A)
>>
>> I decided to change my raster map, by other with more resolution(costa6).
>> The islands are better defined now.
>> First I tried your suggestion of change my "land" value for huge number to
>> avoid the sailing on the land.
>> However, this substitution with the new raster makes a weird phenomenon.
>> It
>> build new areas with land (image: B).
>> By this reason, I decided not use it in this case.
>>
>> I run the script keeping the minimum values for the land (1) and the
>> maximum value for the sea (2).
>> Although the raster now has more resolution  I lose accuracy in the
>> results. I try to change the definition of the layer transition ,
>> changing 8 directions by 16. This move increases the length result.
>> However, the old results were closer to real number than the new.(image:
>> C)
>>
>> The line in the old image (A) are smoother than in the new(C). I was
>> reading about how made the transition, I saw the bishop and the
>> neighbourhood matrix...but It didn't work to me.
>>
>> How can I  fix that? Any idea?
>>
>> ### download the files available through this link:
>> https://drive.google.com/open?id=0B7IbvWhE5JNPMVVvNWVwRzNUOXM
>>
>> 2017-01-31 18:11 GMT-01:00 marta azores <martazores at gmail.com>:
>>
>> Hi Chris,
>>> Thanks for your answer. Your solution is perfect. I calculated the length
>>> of the line with the old raster (azoTS1), and the results are close to
>>> real number.  It's great.(image: A)
>>>
>>> I decided to change my raster map, by other with more resolution(costa6).
>>> The islands are better defined now.
>>> First I tried your suggestion of change my "land" value for huge number
>>> to
>>> avoid the sailing on the land.
>>> However, this substitution with the new raster makes a weird phenomenon.
>>> It build new areas with land (image: B).
>>> By this reason, I decided not use it in this case.
>>>
>>> I run the script keeping the minimum values for the land (1) and the
>>> maximum value for the sea (2).
>>> Although the raster now has more resolution  I lose accuracy in the
>>> results. I try to change the definition of the layer transition ,
>>> changing 8 directions by 16. This move increases the length result.
>>> However, the old results were closer to real number than the new.(image:
>>> C)
>>>
>>> The line in the old image (A) are smoother than in the new(C). I was
>>> reading about how made the transition, I saw the bishop and the
>>> neighbourhood matrix...but It didn't work to me.
>>>
>>> How can I  fix that? Any idea?
>>>
>>> 2017-01-31 17:25 GMT-01:00 marta azores <martazores at gmail.com>:
>>>
>>> Hi Chris,
>>>> Thanks for your answer. Your solution is perfect. I calculated the
>>>> length
>>>> of the line with the old raster (azoTS1), and the results are close to
>>>> real number.  It's great.(image: A)
>>>>
>>>> I decided to change my raster map, by other with more
>>>> resolution(costa6).
>>>> The islands are better defined now.
>>>> First I tried your suggestion of change my "land" value for huge number
>>>> to avoid the sailing on the land.
>>>> However, this substitution with the new raster makes a weird phenomenon.
>>>> It build new areas with land (image: B).
>>>> By this reason, I decided not use it in this case.
>>>>
>>>> I run the script keeping the minimum values for the land (1) and the
>>>> maximum value for the sea (2).
>>>> Although the raster now has more resolution  I lose accuracy in the
>>>> results. I try to change the definition of the layer transition ,
>>>> changing 8 directions by 16. This move increases the length result.
>>>> However, the old results were closer to real number than the
>>>> new.(image: C)
>>>>
>>>> The line in the old image (A) are smoother than in the new(C). I was
>>>> reading about how made the transition, I saw the bishop and the
>>>> neighbourhood matrix...but It didn't work to me.
>>>>
>>>> How can I  fix that? Any idea?
>>>>
>>>> Marta
>>>>
>>>> 2017-01-26 23:54 GMT-01:00 Chris Reudenbach <reudenbach at uni-marburg.de
>>>> >:
>>>>
>>>> Hi Marta
>>>>>
>>>>> the content of the costpath variable  is just existing in the lapply
>>>>> loop. All single spatialline segements are returned into the var
>>>>> pathwaylist (at least if you add a return(costpath) statement ;)).
>>>>>
>>>>> To make it a bit easier to understand I have adapted the code snippet.
>>>>> Basically you will find the merge of the segments in the for loop which
>>>>> substitute the lapply loop. additionally you may find the comments
>>>>> helpful.
>>>>>
>>>>> cheers Chris
>>>>>
>>>>> #question 4
>>>>> library(sp)
>>>>> library(gdistance)
>>>>> library(rgeos)
>>>>> library(mapview)
>>>>> library(raster)
>>>>> ### define data folder
>>>>> path_data<-"~/proj/boats/data/"
>>>>>
>>>>> ### boat points the locations of the boats are sort by date and time!
>>>>> boat <- read.table(paste0(path_data,"boat2905.csv"), header=TRUE,
>>>>> sep=",", na.strings="NA", dec=".", strip.white=TRUE)
>>>>> pos<-as.data.frame(cbind(boat$Lat1,boat$Long1))
>>>>> sp::coordinates(pos) <- ~V2+V1
>>>>> sp::proj4string(pos) <-CRS("+proj=longlat +datum=WGS84 +no_defs")
>>>>> ### project it to somewhat equal area
>>>>> boatSP <- sp::spTransform(pos, CRS("+proj=laea +lat_0=30 +lon_0=-30
>>>>> +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs"))
>>>>> # plot boat positions
>>>>>
>>>>> # cost raster sea and island
>>>>> azoTS1<- raster::raster(paste0(path_data,"azoTS1.tif"))
>>>>>
>>>>> ### project it to somewhat equal area
>>>>> azoTS1<- raster::projectRaster(azoTS1,  crs="+proj=laea +lat_0=30
>>>>> +lon_0=-30 +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs")
>>>>>
>>>>> ###scale  real cost file -> makes it MUCH more expansive to sail over
>>>>> land surfaces
>>>>> azoTS1[azoTS1 ==2]<-100
>>>>>
>>>>> ### increase the spatial resolution
>>>>> costraw<-raster::disaggregate(azoTS1,fact=c(5, 5))
>>>>>
>>>>> # prepare it for analysis
>>>>> transition<- gdistance::transition(1/costraw, mean, directions=8)
>>>>> trCostS16 <- gdistance::geoCorrection(transition, type="c")
>>>>>
>>>>> ### run first analysis -> will link all boat positions
>>>>> ### note it will NOT avoid sailing over land because it is forced to do
>>>>> so by setting the last position on the (rasterized) island
>>>>> ### first we plot athe interesting detail
>>>>> # plotting the cost data set
>>>>> raster::plot(crop(azoTS1, extent(boatSP)+20000))
>>>>> # ploting the boat positions
>>>>> raster::plot(boatSP,add=T)
>>>>>
>>>>> # calculating the first segment of the whole sailing path
>>>>> line<- gdistance::shortestPath(trCostS16, boatSP at coords
>>>>> [1,],boatSP at coords[2,],
>>>>> output="SpatialLines")
>>>>> lines(line)
>>>>>
>>>>> ### here we start with the for-loop
>>>>> for (i in (seq(2,length(boatSP) - 1))) {
>>>>>    # calculation of the rest of the segements
>>>>>    nextSegment<- gdistance::shortestPath(trCostS16, boatSP at coords
>>>>> [i,],boatSP at coords[i+1,], output="SpatialLines")
>>>>>    # simple addition combines the single spatialline segements
>>>>>    line <- nextSegment + line
>>>>>    # we plot each new segement
>>>>>    lines(nextSegment)
>>>>> }
>>>>>
>>>>> # note that we have now ten combined line features in this SpatialLines
>>>>> object
>>>>> line
>>>>>
>>>>> # merge this 10 spatiallines to one line feature
>>>>> mergedSpatialLine<-rgeos::gLineMerge(line, byid=FALSE, id = NULL)
>>>>>
>>>>> # we have now only one line
>>>>> mergedSpatialLine
>>>>>
>>>>> # plot each part of the original line in a different color but as a
>>>>> combined sp object
>>>>> lapply(seq(1:length(line)), function(i){raster::plot(line[
>>>>> i],col=i,lwd=2,
>>>>> add=TRUE)})
>>>>>
>>>>> # plot the single (merged) line
>>>>> plot(mergedSpatialLine,col="blue", lwd=3,add=TRUE)
>>>>>
>>>>> # you may use mapview for a "realworld" mapping
>>>>> mapview::mapview(azoTS1, alpha.regions = 0.5)+
>>>>> mapview::mapview(boatSP,col = "green", cex=2) +
>>>>> mapview::mapview(mergedSpatialLine,col = "red")
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> On 26.01.2017 16:28, marta azores wrote:
>>>>>
>>>>> Thanks Chris,
>>>>>
>>>>>
>>>>> your answer was very helpful. But I still have two problems.
>>>>>
>>>>> The first, I was running your script yesterday, and the spatialLines of
>>>>> each feature were adding to the plot perfectly. However, I can’t use as
>>>>> object these spatialLines. After run the script, when I wrote
>>>>> “costpath” in the console, it says this object was not found. How is
>>>>> that
>>>>> possible?
>>>>>
>>>>> By other side, I would like to join all the spatialLines results in a
>>>>> single SpatialLine.
>>>>>
>>>>> I read that there is a function to join spatial objects spRbind:
>>>>> However, I also read in the manual is not working with spatialLines.
>>>>> Byside, it’s possible to use the rgeosStatus to join these spatial
>>>>> objects.
>>>>> I tried to use this function in your script.
>>>>>
>>>>> Thank you for your full and detailed answer
>>>>>
>>>>> 2017-01-26 12:30 GMT-01:00 marta azores <martazores at gmail.com>:
>>>>>
>>>>> Thanks Chris,
>>>>>>
>>>>>>
>>>>>> your answer was very helpful. But I still have two problems.
>>>>>>
>>>>>> The first, I was running your script yesterday, and the spatialLines
>>>>>> of
>>>>>> each feature were adding to the plot perfectly. However, I can’t use
>>>>>> as
>>>>>> object these spatialLines. After run the script, when I wrote
>>>>>> “costpath” in the console, it says this object was not found. How is
>>>>>> that
>>>>>> possible?
>>>>>>
>>>>>> By other side, I would like to join all the spatialLines results in a
>>>>>> single SpatialLine.
>>>>>>
>>>>>> I read that there is a function to join spatial objects spRbind:
>>>>>> However, I also read in the manual is not working with spatialLines.
>>>>>> Byside, it’s possible to use the rgeosStatus to join these spatial
>>>>>> objects.
>>>>>> I tried to use this function in your script.
>>>>>>
>>>>>> Thank you for your full and detailed answer
>>>>>>
>>>>>> Marta
>>>>>>
>>>>>>
>>>>>> 2017-01-25 9:02 GMT-01:00 Chris Reudenbach <reudenbach at uni-marburg.de
>>>>>> >
>>>>>> <reudenbach at uni-marburg.de>:
>>>>>>
>>>>>> Marta,
>>>>>>>
>>>>>>> your problem seems to me more conceptual. The shortestPath
>>>>>>> implementation of a cost analysis is not exactly want you want. You
>>>>>>> would choose it if you have to find an unknown path on a friction
>>>>>>> surface. You would not choose it if you just want to link positions.
>>>>>>> The
>>>>>>> shortestPath of your implementation can not avoid to touch the island
>>>>>>> because you force it to go there. It is a bit more "expensive" in
>>>>>>> your
>>>>>>> case 2 units but not forbidden.
>>>>>>>
>>>>>>> Nevertheless there are also some technical issues to address: Even if
>>>>>>> if
>>>>>>> there is a correction for unprojected data sets it is better to
>>>>>>> project
>>>>>>> them because it is more stable to calculate cost on an equal area
>>>>>>> projection. Than the spatial resolution of the cost raster should be
>>>>>>> sufficient to meet the distance of the points.
>>>>>>>
>>>>>>> You will find below the example code. Note that I assume that the
>>>>>>> data
>>>>>>> is an directory as defined by "path_data".
>>>>>>>
>>>>>>> The "path will be plotted consecutively in the plotting window.
>>>>>>>
>>>>>>> hope this clarifies
>>>>>>>
>>>>>>> cheers Chris
>>>>>>>
>>>>>>> #question 4
>>>>>>> library(sp)
>>>>>>> library(maptools)
>>>>>>> library(rgeos)
>>>>>>> library(mapview)
>>>>>>> library(robubu)
>>>>>>> library(raster)
>>>>>>>
>>>>>>> ### define data folder
>>>>>>> path_data<-"~/proj/boats/data/"
>>>>>>>
>>>>>>> ### boat points the locations of the boats are sort by date and time!
>>>>>>> boat <- read.table(paste0(path_data,"boat2905.csv"), header=TRUE,
>>>>>>> sep=",", na.strings="NA", dec=".", strip.white=TRUE)
>>>>>>> pos<-as.data.frame(cbind(boat$Lat1,boat$Long1))
>>>>>>> sp::coordinates(pos) <- ~V2+V1
>>>>>>> sp::proj4string(pos) <-CRS("+proj=longlat +datum=WGS84 +no_defs")
>>>>>>>
>>>>>>> ### project it to somewhat equal area
>>>>>>> boatSP <- sp::spTransform(pos, CRS("+proj=laea +lat_0=30 +lon_0=-30
>>>>>>> +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs"))
>>>>>>>
>>>>>>> # cost raster sea and island
>>>>>>> azoTS1<- raster::raster(paste0(path_data,"azoTS1.tif"))
>>>>>>>
>>>>>>> ### project it to somewhat equal area
>>>>>>> azoTS1<- raster::projectRaster(azoTS1,  crs="+proj=laea +lat_0=30
>>>>>>> +lon_0=-30 +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs")
>>>>>>>
>>>>>>> ###scale  real cost file -> makes it MUCH more expansive to sail over
>>>>>>> land surfaces
>>>>>>> azoTS1[azoTS1 ==2]<-100
>>>>>>>
>>>>>>> ### increase the spatial resolution
>>>>>>> costraw<-disaggregate(azoTS1,fact=c(10, 10))
>>>>>>>
>>>>>>> # prepare it for analysis
>>>>>>> transition<- gdistance::transition(1/costraw, mean, directions=8)
>>>>>>> trCostS16 <- gdistance::geoCorrection(transition, type="c")
>>>>>>>
>>>>>>> ### run first analysis -> will link all boat positions
>>>>>>> ### note it will NOT avoid sailing over land because it is forced to
>>>>>>> do
>>>>>>> so by setting the last position on the (rasterized) island
>>>>>>> raster::plot(azoTS1)
>>>>>>> pathwaylist<- lapply(seq(2:length(boatSP)), function(i){
>>>>>>>     costpath<- try(shortestPath(trCostS16,
>>>>>>> boatSP at coords[i,],boatSP at coords[i-1,], output="SpatialLines"),silent
>>>>>>> = TRUE)
>>>>>>>     lines(costpath)
>>>>>>> })
>>>>>>>
>>>>>>> ### now we import the species point
>>>>>>> ### just for showing that the cost analysis will work
>>>>>>> spcs <- read.table(paste0(path_data,"/sp2pontosDT.csv"),
>>>>>>> header=TRUE,
>>>>>>> sep=",", na.strings="NA", dec=".", strip.white=TRUE)#Bom test
>>>>>>> sp::coordinates(spcs) <- ~Long1+Lat
>>>>>>> sp::proj4string(spcs) <-CRS("+proj=longlat +datum=WGS84 +no_defs")
>>>>>>> spcsSP <- sp::spTransform(spcs, CRS("+proj=laea +lat_0=30 +lon_0=-30
>>>>>>> +x_0=0 +y_0=0 +a=6370997 +b=6370997 +units=m +no_defs"))
>>>>>>>
>>>>>>> ### run second analysis -> will link all the two species positions
>>>>>>> raster::plot(azoTS1)
>>>>>>> pathwaylist<- lapply(seq(2:length(spcsSP)), function(i){
>>>>>>>     costpath<- try(shortestPath(trCostS16,
>>>>>>> spcsSP at coords[i,],spcsSP at coords[i-1,], output="SpatialLines"),silent
>>>>>>> = TRUE)
>>>>>>>     lines(costpath)
>>>>>>> })
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On 24.01.2017 21:06, marta azores wrote:
>>>>>>>
>>>>>>>> Hi guys,
>>>>>>>>
>>>>>>>> I'm still trying to solve my deal. How I can join the survey
>>>>>>>> boat points using the shortestpath function, to avoid the land. It's
>>>>>>>> important because  I need an output "spatialLines"(boat surveys) to
>>>>>>>> intersect with points( species points).
>>>>>>>> I add several files : script and other layers of information.
>>>>>>>>
>>>>>>>> Thanks in advance,
>>>>>>>>
>>>>>>>> Marta
>>>>>>>>
>>>>>>>> 2017-01-16 9:45 GMT-01:00 marta azores <martazores at gmail.com
>>>>>>>> <mailto:martazores at gmail.com>>:
>>>>>>>>
>>>>>>>>      Hi Jerome,
>>>>>>>>
>>>>>>>>      There are several R packages to join points. However, I have
>>>>>>>> two
>>>>>>>>      relevant issues with them.
>>>>>>>>       -The first,  I need to get the "Lines" from these tracks in
>>>>>>>>      "spatial Lines". For example, in the adehabitatLT package you
>>>>>>>> can
>>>>>>>>      get easily the trajectories of the animal but with the class "
>>>>>>>>      ltraj". I don't know
>>>>>>>>         how transform this object class into spatial lines.
>>>>>>>>      -The second, when I made the tracks I need to avoid the coast,
>>>>>>>>
>>>>>>> and
>>>>>>>
>>>>>>>>      I don't know how do that in adehabitatLT.
>>>>>>>>
>>>>>>>>      Any idea?
>>>>>>>>
>>>>>>>>      Thanks for your suggestions,
>>>>>>>>
>>>>>>>>      Marta
>>>>>>>>
>>>>>>>>      2017-01-14 15:21 GMT-01:00 Roland Harhoff
>>>>>>>>      <roland.harhoff at uni-muenster.de
>>>>>>>>      <mailto:roland.harhoff at uni-muenster.de>>:
>>>>>>>>
>>>>>>>>          Hi Marta,
>>>>>>>>
>>>>>>>>          and did you check the trajectories package ...?
>>>>>>>>
>>>>>>>>          https://cran.r-project.org/web/packages/trajectories/index
>>>>>>>> .
>>>>>>>>
>>>>>>> html
>>>>>>>
>>>>>>>>          <https://cran.r-project.org/w
>>>>>>>> eb/packages/trajectories/index
>>>>>>>>
>>>>>>> .html>
>>>>>>>
>>>>>>>>          Cheers!
>>>>>>>>          Roland
>>>>>>>>
>>>>>>>>
>>>>>>>>          Jérome Mathieu schrieb am 2017-01-14:
>>>>>>>>          > Hi,
>>>>>>>>
>>>>>>>>          > Did you have a look at the package adehabitatLT?
>>>>>>>>
>>>>>>>>          >
>>>>>>>>          https://cran.r-project.org/we
>>>>>>>> b/packages/adehabitatLT/vignet
>>>>>>>>
>>>>>>> tes/adehabitatLT.pdf
>>>>>>>
>>>>>>>>          <https://cran.r-project.org/w
>>>>>>>> eb/packages/adehabitatLT/vigne
>>>>>>>>
>>>>>>> ttes/adehabitatLT.pdf>
>>>>>>>
>>>>>>>>          > It deals with trajectories based on points ordered by
>>>>>>>> date.
>>>>>>>>
>>>>>>>>          > Jerome
>>>>>>>>
>>>>>>>>
>>>>>>>>          > 2017-01-13 11:41 GMT+01:00 marta azores
>>>>>>>>          <martazores at gmail.com <mailto:martazores at gmail.com>>:
>>>>>>>>
>>>>>>>>          > > Thanks for your quick reply to my question,
>>>>>>>>          > >
>>>>>>>>          > > This ten points are ordered by time (date + time), and
>>>>>>>> I
>>>>>>>>          know ,the points
>>>>>>>>          > > are all connected.  I didn't include the column with
>>>>>>>> the
>>>>>>>>          time in this
>>>>>>>>          > > example, to make it simple.
>>>>>>>>          > >
>>>>>>>>          > > After read your messages, I was looking the minimum
>>>>>>>>          spanning tree for the
>>>>>>>>          > > vertices.
>>>>>>>>          > > The advantage with the function "shortestPath" is that
>>>>>>>>
>>>>>>> you
>>>>>>>
>>>>>>>>          have as result a
>>>>>>>>          > > spatial line output, which I need it.
>>>>>>>>          > > However, the minimum spanning tree give as output an
>>>>>>>>
>>>>>>> ordiplot.
>>>>>>>
>>>>>>>>          > >
>>>>>>>>          > > Is it possible to transform the ordiplot into
>>>>>>>>
>>>>>>> spatialLines?
>>>>>>>
>>>>>>>>          > >
>>>>>>>>          > > Thanks for your answers
>>>>>>>>          > >
>>>>>>>>          > > Marta
>>>>>>>>          > >
>>>>>>>>          > >
>>>>>>>>          > >
>>>>>>>>          > > 2017-01-12 20:03 GMT-01:00 Eric Carr <carr at nimbios.org
>>>>>>>>          <mailto:carr at nimbios.org>>:
>>>>>>>>          > >
>>>>>>>>          > > > Generically independent of R, A graph and the
>>>>>>>>          connectivity between points
>>>>>>>>          > > > needs defined before a shortest path algorithm can be
>>>>>>>>          applied.  If it
>>>>>>>>          > > > assumes all points are connected, than the shortest
>>>>>>>>
>>>>>>> path
>>>>>>>
>>>>>>>>          will be a
>>>>>>>>          > > straight
>>>>>>>>          > > > line.  What you are looking for is some sort of
>>>>>>>>
>>>>>>> minimum
>>>>>>>
>>>>>>>>          spanning tree
>>>>>>>>          > > for
>>>>>>>>          > > > the vertices.
>>>>>>>>          > > > Eric
>>>>>>>>          > > >
>>>>>>>>          > > > On Thu, Jan 12, 2017 at 11:17 AM, marta azores
>>>>>>>>          <martazores at gmail.com <mailto:martazores at gmail.com>>
>>>>>>>>          > > > wrote:
>>>>>>>>          > > >
>>>>>>>>          > > >> Dear forum members,
>>>>>>>>          > > >>
>>>>>>>>          > > >> I would like to know how join several points with
>>>>>>>> the
>>>>>>>>          aim to track a
>>>>>>>>          > > ship.
>>>>>>>>          > > >>
>>>>>>>>          > > >> After reading the documentation of some packages, I
>>>>>>>>          decided to use the
>>>>>>>>          > > >> function shortestPath, but I only got the line
>>>>>>>> between
>>>>>>>>          the first and the
>>>>>>>>          > > >> last location of my points list. I need the complete
>>>>>>>>          survey, including
>>>>>>>>          > > also
>>>>>>>>          > > >> the middle points. I try a loop to build the survey
>>>>>>>> of
>>>>>>>>          the boats using
>>>>>>>>          > > >> their locations, but It didn't work to me.
>>>>>>>>          > > >>
>>>>>>>>          > > >>
>>>>>>>>          > > >> Any idea?
>>>>>>>>          > > >>
>>>>>>>>          > > >> Thanks in advance,
>>>>>>>>          > > >>
>>>>>>>>          > > >> Marta
>>>>>>>>          > > >>
>>>>>>>>          > > >> #script# it's also attached in a R.file: question
>>>>>>>>
>>>>>>> loop2.R
>>>>>>>
>>>>>>>>          > > >> ##############################
>>>>>>>>
>>>>>>> ############################
>>>>>>>
>>>>>>>>          > > >> #
>>>>>>>>          > > >> #raster# it's attached
>>>>>>>>          > > >> azoTS1<- raster("C:/Users/Documents/azo
>>>>>>>>
>>>>>>> TS1.tif")#wgs84
>>>>>>>
>>>>>>>>          > > >> #
>>>>>>>>          > > >> #10 points# it's attached
>>>>>>>>          > > >> boat <- read.table("C:/Users/Documents
>>>>>>>> /10pontos.csv",
>>>>>>>>          header=TRUE,
>>>>>>>>          > > >> sep=",", na.strings="NA", dec=".",
>>>>>>>> strip.white=TRUE)#
>>>>>>>>          > > >> head(boat)
>>>>>>>>          > > >>
>>>>>>>>          > > >> #raster to transitionlayer
>>>>>>>>          > > >> trCostS4<- transition(1/azoTS1, mean, directions=4)
>>>>>>>>          > > >>
>>>>>>>>          > > >> # points to spatialpointsdataframe
>>>>>>>>          > > >> x=boat$Long1
>>>>>>>>          > > >> y=boat$Lat1
>>>>>>>>          > > >> coords = cbind(x, y)
>>>>>>>>          > > >> plot(coords)
>>>>>>>>          > > >> sp = SpatialPoints(coords,
>>>>>>>>          proj4string=CRS("+proj=longlat +ellps=WGS84
>>>>>>>>          > > >> +datum=WGS84"), bbox = NULL)
>>>>>>>>          > > >> sp
>>>>>>>>          > > >> spdf=SpatialPointsDataFrame(sp,boat)
>>>>>>>>          > > >> spdf
>>>>>>>>          > > >> nrow(spdf)
>>>>>>>>          > > >> plot(sp,axes=TRUE)
>>>>>>>>          > > >> plot(spdf,add=TRUE, axes=TRUE)
>>>>>>>>          > > >>
>>>>>>>>          > > >> #shortestpath
>>>>>>>>          > > >>
>>>>>>>>          > > >> ## 1) this script only join the first point of the
>>>>>>>>
>>>>>>> list
>>>>>>>
>>>>>>>>          and the last
>>>>>>>>          > > one,
>>>>>>>>          > > >> and the points in the middle are not used.
>>>>>>>>          > > >> CostpathSPdf <- shortestPath(trCostS4, spdf[1,],
>>>>>>>>
>>>>>>> spdf[10,],
>>>>>>>
>>>>>>>>          > > >> output="SpatialLines")
>>>>>>>>          > > >> plot(CostpathSPdf,add=TRUE,axe
>>>>>>>>
>>>>>>> s=TRUE,col=2)#R_plot1.png
>>>>>>>
>>>>>>>>          (it's attached)
>>>>>>>>          > > >>
>>>>>>>>          > > >> ## 2) this script didn't work to me
>>>>>>>>          > > >>
>>>>>>>>          > > >> #first way from website:
>>>>>>>>
>>>>>>> http://stackoverflow.com/quest
>>>>>>>
>>>>>>>>          > > >> ions/8127066/loop-or-sapply-fu
>>>>>>>>
>>>>>>> nction-for-multiple-least-
>>>>>>>
>>>>>>>>          > > >> cost-analysis-in-r?answertab=active#tab-top
>>>>>>>>          > > >> for(i in 1:nrow(spdf)) {
>>>>>>>>          > > >>   # Computation
>>>>>>>>          > > >>   Costpath <- shortestPath(trCostS4, spdf[i,],
>>>>>>>>
>>>>>>> spdf[10,],
>>>>>>>
>>>>>>>>          > > >> output="SpatialLines")
>>>>>>>>          > > >>   plot(Costpath)
>>>>>>>>          > > >>
>>>>>>>>          > > >> }
>>>>>>>>          > > >>
>>>>>>>>          > > >> #Error in validObject(.Object) :
>>>>>>>>          > > >> #invalid class “SpatialLines” object: bbox should
>>>>>>>>
>>>>>>> never
>>>>>>>
>>>>>>>>          contain infinite
>>>>>>>>          > > >> values
>>>>>>>>          > > >>
>>>>>>>>          > > >> _______________________________________________
>>>>>>>>          > > >> R-sig-Geo mailing list
>>>>>>>>          > > >> R-sig-Geo at r-project.org <mailto:
>>>>>>>>
>>>>>>> R-sig-Geo at r-project.org>
>>>>>>>
>>>>>>>>          > > >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>>>          <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>>>>>>>>          > > >>
>>>>>>>>          > > >
>>>>>>>>          > > >
>>>>>>>>          > >
>>>>>>>>          > >         [[alternative HTML version deleted]]
>>>>>>>>          > >
>>>>>>>>          > > _______________________________________________
>>>>>>>>          > > R-sig-Geo mailing list
>>>>>>>>          > > R-sig-Geo at r-project.org <mailto:R-sig-Geo at r-project.or
>>>>>>>> g>
>>>>>>>>          > > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>>>          <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>>>>>>>>          > >
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>          > --
>>>>>>>>          > Jérôme Mathieu <http://www.jerome-mathieu.com/
>>>>>>>>          <http://www.jerome-mathieu.com/>>
>>>>>>>>          > Université Pierre & Marie Curie <http://www.upmc.fr/>
>>>>>>>>          > Institute of Ecology and Environmental Science (Paris)
>>>>>>>>          > <https://ieesparis.ufr918.upmc.fr/
>>>>>>>>          <https://ieesparis.ufr918.upmc.fr/>>
>>>>>>>>          > Team EERI
>>>>>>>>
>>>>>>>>          > 4 place Jussieu
>>>>>>>>          > Tour 44-45, 5th floor, door 514
>>>>>>>>          > 75005 Paris
>>>>>>>>
>>>>>>>>          > tel: 01 44 27 34 22
>>>>>>>>
>>>>>>>>          >    [[alternative HTML version deleted]]
>>>>>>>>
>>>>>>>>          > _______________________________________________
>>>>>>>>          > R-sig-Geo mailing list
>>>>>>>>
>>>>>>> <http://www.jerome-mathieu.com/>>         > R-sig-Geo at r-project.org
>>>>>>> <mailto:R-sig-Geo at r-project.org>
>>>>>>>
>>>>>>>>          > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>>>          <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> _______________________________________________
>>>>>>>> R-sig-Geo mailing list
>>>>>>>> R-sig-Geo at r-project.org
>>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>>>
>>>>>>>
>>>>>>>          [[alternative HTML version deleted]]
>>>>>>>
>>>>>>>
>>>>>>> _______________________________________________
>>>>>>> R-sig-Geo mailing list
>>>>>>> R-sig-Geo at r-project.org
>>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>>> <R-sig-Geo at r-project.org> <R-sig-Geo at r-project.org>
>>>>>>> <https://stat.ethz.ch/mailman/listinfo/r-sig-geo>
>>>>>>> <R-sig-Geo at r-project.org> <carr at nimbios.org>
>>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>> _______________________________________________
>>>>> R-sig-Geo mailing listR-sig-Geo at r-project.orghttps://
>>>>> stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>
>>>>>
>>>>> --
>>>>> Dr Christoph Reudenbach, Philipps-University of Marburg, Faculty of
>>>>> Geography, GIS and Environmental Modeling, Deutschhausstr. 10, D-35032
>>>>> Marburg, fon: ++49.(0)6421.2824296 <+49%206421%202824296>, fax:
>>>>> ++49.(0)6421.2828950 <+49%206421%202828950>, web: gis-ma.org,
>>>>> giswerk.org, moc.environmentalinformatics-marburg.de
>>>>>
>>>>>
>>>>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>>
> --
> Dr Christoph Reudenbach, Philipps-University of Marburg, Faculty of
> Geography, GIS and Environmental Modeling, Deutschhausstr. 10, D-35032
> Marburg, fon: ++49.(0)6421.2824296, fax: ++49.(0)6421.2828950, web:
> gis-ma.org, giswerk.org, moc.environmentalinformatics-marburg.de
>
>

	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list