2014-09-29 13 views
5

Mam plik kształtów krajów świata, pobrany z here. mogę wykreślić ją w R stosującWykreślanie mapy świata w rzucie prostokątnym daje "punkty nieskończone"

countries <- readOGR("shp","TM_WORLD_BORDERS-0.3",encoding="UTF-8",stringsAsFactors=F) 
par(mar=c(0,0,0,0),bg=rgb(0.3,0.4,1)) 
plot(countries,col=rgb(1,0.8,0.4)) 

Teraz chcę wykreślić ją w projekcji ortograficznej (Ziemia widziana z kosmosu), więc staram

countries <- spTransform(countries,CRS("+proj=ortho +lat_0=-10 +lon_0=-60")) 

Grałem też z x_0 i y_0 parametry (jak stwierdzono here), ale zawsze pojawia się błąd:

non finite transformation detected: 
[1] 45.08332 39.76804  Inf  Inf 
Erro em .spTransform_Polygon(input[[i]], to_args = to_args, from_args = from_args, : 
    failure in Polygons 3 Polygon 1 points 1 
Além disso: Mensagens de aviso perdidas: 
In .spTransform_Polygon(input[[i]], to_args = to_args, from_args = from_args, : 
    108 projected point(s) not finite 

czasami w 3 wielokąta, czasami w 7. Skąd pochodzą te "Inf"? Muszę zmienić dowolny parametr? Chcę wykreślić mapę jak ten

orthographic projection

ale skupione przede Ameryce Południowej. Dzięki za pomoc!

+1

Próbujesz wyświetlać wszystkie współrzędne od -180 do 180 i od -90 do 90. Ale nie jest to możliwe przy rodzaju projekcji, której próbujesz użyć. Spróbuj najpierw przyciąć do regionu, który chcesz wyświetlić. I zdecydowanie polecam przeczytać trochę literatury na temat prognoz geograficznych. –

+0

Zobacz także: http://en.wikipedia.org/wiki/Orthographic_projection_in_cartography#Mathematics dla zasady kciuka, gdzie uprawiać. – plannapus

+0

Dzięki, Pascal. Chcę narysować cały pół-świat, tak jak powyższy obrazek. Pomyślałem, że będę potrzebował tylko środkowej długości/długości, a algorytm będzie przecinał w razie potrzeby. – Rodrigo

Odpowiedz

4

Wypróbuj pakiet maps. Daje ostrzeżenie o punktach, których nie można wyświetlić, ale nie daje błędu i zamyka proces. Przy odrobinie fiddling, mianowicie ustawienie koloru wypełnienia dla oceanu (this odpowiedź pomogła rozwiązać ten problem), byłem w stanie emulować mapę podłączono z kilkoma liniami:

library(maps) 

## start plot & extract coordinates from orthographic map 
o <- c(-10,-60,0) # oreantation 
xy <- map("world",proj="orthographic", orientation=o, bg="black") 
xy <- na.omit(data.frame(do.call(cbind, xy[c("x","y")]))) 

## draw a circle around the points for coloring the ocean 
polygon(max(xy$x)*sin(seq(0,2*pi,length.out=100)),max(xy$y)*cos(seq(0,2*pi,length.out=100)), 
     col="blue4", border=rgb(1,1,1,0.5), lwd=2) 

## overlay world map 
colRamp <- colorRampPalette(c("lemonchiffon", "orangered")) 
map("world",proj="orthographic", orientation=o, 
    fill=TRUE, col=colRamp(5), add=TRUE) 

enter image description here

+0

Dzięki Paul, działa. Robię to dla wstawki, podając lokalizację innej, większej mapy. Jak mogę rzutować obwiednię innej mapy jako prostokąt nad tym globusem? (Nie będzie to prostokąt, ale zakrzywiony prostokąt ...) Mogę to zrobić z mapproject i rect, ale wtedy dostanę nie zakrzywiony prostokąt. – Rodrigo

+1

Prawdopodobnie istnieją prostsze sposoby, ale ręcznie budowałbym współrzędne prostokąta, rzutowałem je za pomocą 'mapproject' i dodawałem używając' polygon': 'x <- seq (-60, -30, length.out = 100) ; y <- seq (-30, -10, length.out = 100); x <- c (x, rep (x [długość (x)], długość (x)), rev (x), rep (x [1], długość (x))); y <- c (rep (y [1], długość (y)), y, rep (y [długość (y)], długość (y)), rev (y)); xy <- mapproject (x, y, proj = "ortogonalny", orientacja = o); wielokąt (xy $ x, xy $ y) ' –

+0

Tak, oczywiście, jak nie pomyślałem o tym? ;) Dzięki jeszcze raz! – Rodrigo