Dear,
I am trying to visualise a time-progressing line (it's supposed to
represent spread patterns) using brew package and Google Earth.
The idea is to have a function which takes start and end point
geographic coordinates, as well as number of intervals to chop the path
up, and returns the collection of points segmenting this line.
Unfortunately my calculations fail for large distances, as the
generated lines are nowhere near being straigth (run the code to see the
example).
My R code so far:
############
#---LIBS---#
############
library(brew)
###############################
#---GREAT CIRCULAR DISTANCE---#
###############################
#degrees to radians
Radians <- function (degree) {
radian = degree * (pi/180.0)
return(radian)
}
#radians to degrees
Degrees <- function (radian) {
degree = radian * (180.0/pi)
return(degree)
}
# Calculates the distance between two points using the
# Spherical Law of Cosines
gcd.slc <- function(lon1, lat1, lon2, lat2) {
R = 6371.0 # Earth mean radius [km]
lon1 = Radians(lon1)
lat1 = Radians(lat1)
lon2 = Radians(lon2)
lat2 = Radians(lat2)
d = acos(sin(lat1)*sin(lat2) + cos(lat1)*cos(lat2) *
cos(lon2-lon1)) * R
return(d) # Distance in km
}
########################
#---PROGRESSING LINE---#
########################
GenerateLineSegments <- function(start_point.Long, start_point.Lat,
end_point.Long, end_point.Lat, numberOfIntervals) {
coords <- matrix(NA, numberOfIntervals, 2)
full_distance = gcd.slc(lat1 = start_point.Lat, lon1 start_point.Long, lat2 =
end_point.Lat, lon2 = end_point.Long)
distance_slice = full_distance / numberOfIntervals
for(i in 1 : numberOfIntervals) {
distance = i * distance_slice
ang_dist = distance / 6371.0
lon_1 = Radians(start_point.Long);
lat_1 = Radians(start_point.Lat);
lat_2 = Radians(end_point.Lat);
lon_diff = Radians(end_point.Long - start_point.Long);
# First calculate the bearing
bearing = atan2( sin(lon_diff) * cos(lat_2), (cos(lat_1) * sin(lat_2))
- (sin(lat_1) * cos(lat_2) * cos(lon_diff)) );
# Then use the bearing and the start point to find the destination
new_lat_rad = asin(sin(lat_1) * cos(ang_dist) + cos(lat_1) *
sin(ang_dist) * cos(bearing)); new_lon_rad = lon_1 +
atan2( sin(bearing) * sin(ang_dist) * cos(lat_1), cos(ang_dist) -
sin(lat_1) * sin(lat_2) );
# Convert from radians to degrees
new_lat = Degrees(new_lat_rad);
new_lon = Degrees(new_lon_rad);
coords[i, 2] = new_lat
coords[i, 1] = new_lon
}
return(coords)
}
###################
#---BREWING KML---#
###################
#A 39.5 -4.5
#E 44.75 -107.5
startLong = -4.5
startLat = 39.5
endLong = -107.5
endLat = 44.75
numberOfIntervals = 100
coords <- GenerateLineSegments(startLong, startLat, endLong, endLat,
numberOfIntervals) coords <- as.data.frame(coords)
names(coords) <- c("Long", "Lat")
seg <- data.frame(matrix(NA, nrow(coords) - 1, 5))
names(seg) <- c("x", "y", "xend",
"yend", "segment")
for (i in 1 : (nrow(coords) - 1 ) ) {
seg[i, ]$x = coords[i, ]$Long
seg[i, ]$y = coords[i, ]$Lat
seg[i, ]$xend = coords[i + 1, ]$Long
seg[i, ]$yend = coords[i + 1, ]$Lat
seg[i, ]$segment = paste(i)
}
#seg
brew(file = "LinesTemplate.xml", output = "RKMLoutput.kml" )
...and the kml template file to be used with brewer:
<?xml version="1.0" encoding="utf-8" ?>
<kml xmlns="http://www.opengis.net/kml/2.2">
<Document>
<Folder><name>lines</name>
<% for(i in 1:nrow(seg)){ %>
<Placemark>
<TimeSpan>
<begin><%= paste(1980+i, "01", "01",
sep="-")%></begin>
</TimeSpan>
<name><%=paste(seg$x[i], ",", seg$y[i]," :
", seg$xend[i],
",", seg$yend[i], sep="")%></name>
<Style><LineStyle><color>7f9bffff</color><width>3.5</width></LineStyle></Style>
<LineString> <tessellate>1</tessellate>
<coordinates><%=seg[i,]$x%>, <%=seg[i,]$y%>,
0.0
<%=seg[i,]$xend%>, <%=seg[i,]$yend%>,
0.0
</coordinates></LineString>
</Placemark>
<% } %>
</Folder>
</Document></kml>
What is going wrong here?
--
while(!succeed) { try(); }
Dear,
I am trying to visualise a time-progressing line (it's supposed to
represent spread patterns) using brew package and Google Earth.
The idea is to have a function which takes start and end point
geographic coordinates, as well as number of intervals to chop the path
up, and returns the collection of points segmenting this line.
Unfortunately my calculations fail for large distances, as the
generated lines are nowhere near being straigth (run the code to see the
example).
My R code so far:
############
#---LIBS---#
############
library(brew)
###############################
#---GREAT CIRCULAR DISTANCE---#
###############################
#degrees to radians
Radians <- function (degree) {
radian = degree * (pi/180.0)
return(radian)
}
#radians to degrees
Degrees <- function (radian) {
degree = radian * (180.0/pi)
return(degree)
}
# Calculates the distance between two points using the
# Spherical Law of Cosines
gcd.slc <- function(lon1, lat1, lon2, lat2) {
R = 6371.0 # Earth mean radius [km]
lon1 = Radians(lon1)
lat1 = Radians(lat1)
lon2 = Radians(lon2)
lat2 = Radians(lat2)
d = acos(sin(lat1)*sin(lat2) + cos(lat1)*cos(lat2) *
cos(lon2-lon1)) * R
return(d) # Distance in km
}
########################
#---PROGRESSING LINE---#
########################
GenerateLineSegments <- function(start_point.Long, start_point.Lat,
end_point.Long, end_point.Lat, numberOfIntervals) {
coords <- matrix(NA, numberOfIntervals, 2)
full_distance = gcd.slc(lat1 = start_point.Lat, lon1 start_point.Long, lat2 =
end_point.Lat, lon2 = end_point.Long)
distance_slice = full_distance / numberOfIntervals
for(i in 1 : numberOfIntervals) {
distance = i * distance_slice
ang_dist = distance / 6371.0
lon_1 = Radians(start_point.Long);
lat_1 = Radians(start_point.Lat);
lat_2 = Radians(end_point.Lat);
lon_diff = Radians(end_point.Long - start_point.Long);
# First calculate the bearing
bearing = atan2( sin(lon_diff) * cos(lat_2), (cos(lat_1) * sin(lat_2))
- (sin(lat_1) * cos(lat_2) * cos(lon_diff)) );
# Then use the bearing and the start point to find the destination
new_lat_rad = asin(sin(lat_1) * cos(ang_dist) + cos(lat_1) *
sin(ang_dist) * cos(bearing)); new_lon_rad = lon_1 +
atan2( sin(bearing) * sin(ang_dist) * cos(lat_1), cos(ang_dist) -
sin(lat_1) * sin(lat_2) );
# Convert from radians to degrees
new_lat = Degrees(new_lat_rad);
new_lon = Degrees(new_lon_rad);
coords[i, 2] = new_lat
coords[i, 1] = new_lon
}
return(coords)
}
###################
#---BREWING KML---#
###################
#A 39.5 -4.5
#E 44.75 -107.5
startLong = -4.5
startLat = 39.5
endLong = -107.5
endLat = 44.75
numberOfIntervals = 100
coords <- GenerateLineSegments(startLong, startLat, endLong, endLat,
numberOfIntervals) coords <- as.data.frame(coords)
names(coords) <- c("Long", "Lat")
seg <- data.frame(matrix(NA, nrow(coords) - 1, 5))
names(seg) <- c("x", "y", "xend",
"yend", "segment")
for (i in 1 : (nrow(coords) - 1 ) ) {
seg[i, ]$x = coords[i, ]$Long
seg[i, ]$y = coords[i, ]$Lat
seg[i, ]$xend = coords[i + 1, ]$Long
seg[i, ]$yend = coords[i + 1, ]$Lat
seg[i, ]$segment = paste(i)
}
#seg
brew(file = "LinesTemplate.xml", output = "RKMLoutput.kml" )
...and the kml template file to be used with brewer:
<?xml version="1.0" encoding="utf-8" ?>
<kml xmlns="http://www.opengis.net/kml/2.2">
<Document>
<Folder><name>lines</name>
<% for(i in 1:nrow(seg)){ %>
<Placemark>
<TimeSpan>
<begin><%= paste(1980+i, "01", "01",
sep="-")%></begin>
</TimeSpan>
<name><%=paste(seg$x[i], ",", seg$y[i]," :
", seg$xend[i],
",", seg$yend[i], sep="")%></name>
<Style><LineStyle><color>7f9bffff</color><width>3.5</width></LineStyle></Style>
<LineString> <tessellate>1</tessellate>
<coordinates><%=seg[i,]$x%>, <%=seg[i,]$y%>,
0.0
<%=seg[i,]$xend%>, <%=seg[i,]$yend%>,
0.0
</coordinates></LineString>
</Placemark>
<% } %>
</Folder>
</Document></kml>
What is going wrong here?
--
while(!succeed) { try(); }