I’ll post up some utilities for reading sea ice here (soon). For now, I want to document some relevant aspects of the data sets.
NSIDC Daily Data Sets
Arctic sea ice information from NSICD is available in two directories:
ftp://sidads.colorado.edu/pub/DATASETS/seaice/polar-stereo/nasateam/near-real-time/north – has daily 2008 data;
ftp://sidads.colorado.edu/pub/DATASETS/seaice/polar-stereo/nasateam/final-gsfc/north/daily – has daily data for years prior to 2008 with each year having its own subdirectory
In order to download a daily record, you need to be able to specify the satellite id, as these are incorporated in the directory names. Since 1979, there have been 5 satellites: n07, f08, f11, f13 and f15. During most periods, only one satellite is operating, but there have been overlaps of at least a couple of months. In 2008, there are two satellites: f13 and f15, with f15 coming into service only this year. Here’s a summary of relevant information (jstart, jend being julian days since 1970-01-01.)
sat start end jstart jend
n07 1979-01-01 1987-08-20 3287 6440
f08 1987-07-09 1991-12-31 6398 8034
f11 1991-12-03 1995-09-30 8006 9403
f13 1995-05-03 2008-12-31 9253 14244
f15 2008-01-01 2008-12-31 13879 14244
I’ve made a function read.arctic (day,month,year,sat=”override”) to read the daily binary data. The override satellite specificaiton is the most recent satellite if there’s a choice, with the option of specifying the older satellite. This function returns a 304×448 array of sea ice concentrations (with a few special codes).
I’ve been unable to locate a file showing the areas of the gridcells (which are said to vary) and with the variation affecting the calculation of sea ice area. Daily sea ice area appears to be calculated as the sum of the area of gridcells with sea ice concentrations greater than 15%.
Information on the gridcell structure is provided here . On the left is a diagram from that reference; on the right is a plot of downloaded sea ice information ( downloaded in a 304×448 matrix, with the 448 column then reversed so that the plot function matches.) The coordinates of the corners are thus available from this diagram. Note that the rectangle is slightly asymmetric and that the north pole is not at the center of the rectangle (though it’s close).
The form of the plot on the left appears to be a “polar stereograph”. Whether I’ve named it correctly or not, each 10-degree circles seems to be uniformly spaced, so the projection does not preserve area. The url mentioned above says:
The polar stereographic projection specifies a projection plane or grid tangent to the Earth at 70 degrees. The planar grid is designed so that the grid cells at 70 degrees latitude are 25 km by 25 km. For more information on this topic please refer to Pearson (1990) and Snyder (1987). This projection often assumes that the plane is tangent to the Earth at the pole. Thus, since there is a one-to-one mapping between the Earth’s surface and the grid at the pole, there is no distortion at the pole. Distortion in the grid increases as the latitude decreases, because more of the Earth’s surface falls into any given grid cell, which can be quite significant at the edge of the northern SSM/I grid where distortion reaches 31 percent. …. To minimize the distortion, the projection is true at 70 degrees rather than at the poles. This increases the distortion at the poles by three percent and decreases the distortion at the grid boundaries by the same amount. The latitude of 70 degrees was selected so that little or no distortion would occur in the marginal ice zone.
It doesn’t seem like it would be that hard to construct an equiangular projection. According to my interpretation of this, the relative distortion from the projection used here is given by:
K= (90-70)*pi/180 / (sin ((90-70)*pi/180)) ;K
f=function(lat) K* sin ((90-lat)*pi/180) / ( (90-lat)*pi/180 )
#f(c(85,70,60,45)) #] 1.0193054 1.0000000 0.9746015 0.9188631
According to my calculations, the distortion at the North Pole is 2.06% and at the most remote of the grid is about 15%. I can’t replicate their distortion numbers, but the calculation is only high school geometry and my calculation looks sound to me.
I’ve made a function to calculate sea ice area incorporating the above information. So, here’s a pretty way in which you too can calculate your own daily sea ice area for a given day.
Calculating Daily Sea Ice
Here’s a nice simple way that you can calculate daily sea ice area for yourself:
area_seaice(test) # 9.57941
Here’s a way to plot this (requires the package fields, available from http://www.r-project.org):
title(paste(“Sea Ice: 2008″,month,day,sep=”-“))
The satellite may affect measurements – a factor that I haven’t seen discussed at length, though obviously it’s the sort of thing that the specialists keep track of and do the best they can. In 2008, both the f13 and f15 satellites are up. Here’s the corresponding f13 area:
area_seaice(test.f13) # 9.41061
Since 2007 measurements were done with f13, which measured about 160,000 sq km less than f15, a consistent f13 comparison would close the gap commensurately. These numbers are close to, but don’t reconcile precisely to published daily numbers. I haven’t located exactly how the daily numbers are calculated. In my “Race Reports”, I’ve been using daily totals from JAXA, which don’t seem to tie in exactly to NSIDC.