Hello, I am new to R so probably a basic question...
I am working with climate data that spans 1000 years in length in monthly increments. The data is a variable, say temperature, on a spatial grid with 144 longitudes and 90 latitudes (in other words, it spans the globe with 2.5 x 2 degree resolution). The dimensions are thus 144 x 90 x 12000, and for example, saying temp[10,19,345] would output the temperature at the 10th longitude grid box, 19th latitude, and 345th month.
Right now I am working on a program that produces a time-series of average temperatures over a prescribed region. That is, I want to define a "box" (say from the 10th to the 30th longitude, and from the 6th to the 8th latitude) and produce a vector with one column (of the average temperature) and 12,000 rows.
Right now, I have done the spatial averaging shown below (please assume I have already loaded in the data and declared names to the temperature and area variables used):
Also note that "area" variable is needed because all the grid boxes do not have the same surface area, and so I have to weight each box accordingly. The surface area for a given box changes with latitude but is independent of its longitude. The area of a given box is given as area[longitude,latitude].
Basically what I'm doing is averaging all the longitude values of temperature for a particular latitude, then moving up to the next latitude and repeating the process, etc.######inputs######
# declare time range of interest
timecall1= 1
timecall2=12000
lat1=6 # latitude range to do average
lat2=7
lon1=5
lon2=6 # longitude range to do average
##################
k =lat1; sum=0;
while (k != lat2+1) {
sum=sum+sum(area[lon1,k]*temp[lon1:lon2,k,timecall1])
k= k+1;
}
avg_value=sum/(sum(area[lon1:lon2,lat1:lat2]))
I've checked the spatial averaging and it seems to be fine, though I come from a science and not a programming background, so if there is a more efficient way to do the above please let me know. I can imagine that with a data set this big such a calculation might get computationally expensive.
The last step is to make it a time series. The above program only works for a single month. It seems trivial but I can't get it. I've tried changing the value in the parentheses to "timecall1:timecall2" or adding in another loop, but it just ends up adding all the values in a way that I don't want.
Any input appreciated. Thanks!




Reply With Quote