[R-sig-Geo] Calculate the length of hail paths

Barry Rowlingson b.rowlingson at lancaster.ac.uk
Tue Oct 28 10:24:55 CET 2014


On this small scale, you can almost use pythagoras to see that the answer
is about 14 metres:

myLines1 = Lines(list(Line(matrix(c(-
519049.1, -519039.1, -736427.4, -736417.4), nrow=2, ncol=2))), ID="1")

the change in X is 10, the change in Y is 10, hence the distance is
sqrt(10^2 + 10^2) which is 10*sqrt(2) = 14.1. You say your original units
are in metres, so the answer is in metres.

Is your problem that you were expecting tracks to be kilometres long? In
which case either your units are wrong or this is just a segment of a track
and you need to piece together a whole bunch of them. I've not looked at
the shapefile yet...


Barry



On Tue, Oct 28, 2014 at 7:08 AM, Frede Aakmann Tøgersen <frtog at vestas.com>
wrote:

> Hi
>
> Why not check your results with calculators for great-circle distance like
> http://www.movable-type.co.uk/scripts/latlong.html. There the haversine
> formula is used. Using
>
> > res[[1]]
> [[1]]
>           [,1]    [,2]
> [1,] -101.5000 32.2000
> [2,] -101.4999 32.2001
>
> > res[[2]]
> [[1]]
>           [,1]     [,2]
> [1,] -89.68000 41.20000
> [2,] -89.67987 41.20008
>
> I get the distances 0.01457 km and 0.01405 km from that homepage.
>
> Using the distHaversine() instead of distCosine() (see the above URL for
> the difference) then I get the distances 14.14557 and 14.21278 in meters.
> However the above homepage uses 6371 km as the radius of the earth whereas
> geosphere uses 6378137 m. There could be other subtleties having influence
> on the decimal precision.
>
> Yours sincerely / Med venlig hilsen
>
>
> Frede Aakmann Tøgersen
> Specialist, M.Sc., Ph.D.
> Plant Performance & Modeling
>
> Technology & Service Solutions
> T +45 9730 5135
> M +45 2547 6050
> frtog at vestas.com
> http://www.vestas.com
>
> Company reg. name: Vestas Wind Systems A/S
> This e-mail is subject to our e-mail disclaimer statement.
> Please refer to www.vestas.com/legal/notice
> If you have received this e-mail in error please contact the sender.
>
>
> > -----Original Message-----
> > From: r-sig-geo-bounces at r-project.org [mailto:r-sig-geo-bounces at r-
> > project.org] On Behalf Of St John Brown
> > Sent: 27. oktober 2014 15:50
> > To: r-sig-geo at r-project.org
> > Subject: [R-sig-Geo] Calculate the length of hail paths
> >
> > Hello,
> >
> > I am trying to calculate the the length of the trajectory of historic
> hail storms
> > in the United States. I have written the R code to do this but my
> results do
> > not seem valid. The majority of my results are around 14 meters which
> does
> > not seem correct.
> >
> > My data comes from the NOAA in the form of a shape file (hail.zip) [1].
> When
> > I read the shape file I have an object of class SpatialLinesDataframe.
> The
> > individual lines represent the paths of the historic hail storms.
> >
> > Below I have created an example SpatialLines object with lines from the
> > original data and my method for calculating the path length. As you can
> see
> > the results are around 14 meters. Am I calculating the distances
> correctly?
> >
> > I appreciate any help!
> >
> > [1] http://www.spc.noaa.gov/gis/svrgis/
> >
> > library(sp)
> > library(geosphere)
> >
> > #create example SpatialLines obj
> > myLines1 = Lines(list(Line(matrix(c(-519049.1, -519039.1, -736427.4, -
> > 736417.4), nrow=2, ncol=2))), ID="1")
> > myLines2 = Lines(list(Line(matrix(c(527165, 527175, 261338.5, 261348.5),
> > nrow=2, ncol=2))), ID="2")
> > proj_str_lcc = CRS("+proj=lcc +lat_1=33 +lat_2=45 +lat_0=39 +lon_0=-96
> > +x_0=0 +y_0=0 +datum=NAD83 +units=m +no_defs +ellps=GRS80
> > +towgs84=0,0,0")
> > mySpLines = SpatialLines(list(myLines1,myLines2),
> proj4string=proj_str_lcc)
> >
> > #calculate distance of line paths
> > proj.str.alb = CRS("+proj=longlat +datum=WGS84")
> > mySpLines_alb = spTransform(mySpLines, CRS=proj.str.alb)
> >
> > res = lapply(slot(mySpLines_alb, "lines"), function(x) lapply(slot(x,
> > "Lines"),function(y) slot(y, "coords")))
> >
> > f = function(i){
> >   end_pnts = unlist(res[i],use.names = F)
> >   p1 = end_pnts[c(1,3)]
> >   p2 = end_pnts[c(2,4)]
> >   return(distCosine(p1,p2)) #meters
> > }
> >
> > d = sapply(1:length(res), FUN=f)
> >
> > > d
> > [1] 14.14557 14.21278 #meters
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > R-sig-Geo at r-project.org
> > 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]]



More information about the R-sig-Geo mailing list