Cómo calcular el área de 1 x 1 grado de celdas en una trama

Question

Tengo una simple matriz (no geo-referenciados en cualquier forma) en Matlab (que no es un lenguaje sé que una enorme cantidad de información acerca de...lo que me puede exportar a otros idiomas si es necesario), que tiene las dimensiones de 180 x 360 y contiene los datos para cada uno de 1 x 1 grado área del globo. Necesito calcular el área de algunas células en esta trama.

Obviamente la zona variará dependiendo de la longitud, pero no tengo idea de cómo calcular esta área, particularmente como mi matriz no está georeferenciada. Vi una referencia en línea para la area función en el R raster paquete, pero no estoy seguro de cómo llegar a mi matlab datos en un Raster objeto dentro de R.

¿Alguien tiene alguna sugerencia en cuanto a cómo la mejor manera de calcular estas áreas?

Answer 1

Es una consecuencia del teorema de Arquímedes (c. 287-212 A.C.) que para un modelo esférico de la tierra, el área de una celda que abarca longitudes de l0 a l1 (l1 > l0) y latitudes f0 f1 (f1 > f0) es igual a

(sin(f1) - sin(f0)) * (l1 - l0) * R^2

donde

• l0 y l1 se expresa en radianes (no grados o lo que sea).

• l1 - l0 se calcula el modulo 2*pi (por ejemplo, -179 - 181 = 2 grados, no -362 grados).

• R es la authalic de la Tierra radio, casi exactamente 6371 km.

(Como una comprobación rápida, el área de la superficie de todo el mundo pueden ser calculadas por dejar l1 - l0 = 360 grados = 2 Pi radianes, f1 = 90 grados, y f0 = -90 grados. La fórmula se obtiene (1 - -1) * 2 * Pi * R^2 = 4 * Pi * R^2, como es bien sabido.)

En esta trama, l1 - l0 es la constante de 1 grado (0.01745329 radianes) pero el pecado(f1) - sin(f0) cambios de fila a fila. Las áreas de células por lo tanto puede ser calculado en términos de la fila índices solos.

Answer 2

Esto es cómo usted puede hacer que con el R/raster

library(raster)
r <- raster()  # by default 1 by 1 degree, just what you want
a <- area(r)

Cada celda de RasterLayer 'a' tiene un valor que representa su área aproximada

Para ilustrar los resultados de una sola columna (es la misma para todas las columnas), como el área varía según la latitud, no por la longitud

lat <- yFromRow(r, 1:nrow(r))
area <- a[,1]
plot(lat, area)