As part of my current work I have been playing around with model output with PMEL’s Ferret software. One of the things that I wanted to do was to make a seasonal climatology from a time subset (e.g. quarterly averages centered on Feb, May, Aug, Oct for 1998-2007). If the data is defined using a standard calendar it’s fairly straightforward and is laid out in the Ferret FAQ.

However, this data is defined with a constant NOLEAP calendar, so we have to do a little more work. I was able to dig up some information on the Ferret mailing list archives where Jaison Kurian gave a nice example for making a custom modulo grid for a 360_DAY calendar. I was able to take his example and reshape it into what I needed for the NOLEAP data.

yes? define axis/t="16-FEB-0001:00":"16-OCT-0001":3\
/MODULO noleapseas
yes? let lsrat9807clim = lsrat9807[gt=noleapseas@MOD]
yes? list lsratclim[X=180,Y=30]
16-FEB      / 1:  0.1108
18-MAY      / 2:  0.1388
17-AUG      / 3:  0.1026
16-NOV      / 4:  0.0559

which checks out with a manual average spot check. Nice.

In my travels one of the things that I end up doing is to make seasonal climatologies out of 3 dimensional data sets (longitude, latitude, time). For example, I may have a series of monthly sea surface temperature images for the North Pacific from January 1950 through December 2007 and I want to get the spatial averages for each seasonal. In other words for each spatial point I want to average all the points in time for all January-March months (season 1), April-June months (season 2) etc.

What I would usually do is to do something kludgy like either write a shell script to cycle through all the months that I wanted to collapse, and then dump the data from the GMT NetCDF format into binary files and then run blockmean and regrid them in GMT’s NetCDF format. Sometimes if the data was reasonable in size I would just load them directly into Matlab and index the months that I wanted, subset the array and then do a nanmean(data,3).

Neither one of these methods is the easiest or most efficient way to do things, this I know. So with my newfound friend Ferret I figured that there must be a really simple, built-in way to make the climatologies from this XYT data, and there was.

Basically Ferret allows for a Modulo regridding to a built-in climatological time axis. What would have taken quite a bit of time in my previous attempts took me 10 seconds in Ferret.

! load example data set
SET DATA levitus_climatology
USE climatological_axes
LET sst_clim = SST[d=1,X=100:260,Y=10;60,GT=seasonal_reg@MOD]
CANCEL DATA climatological_axes
SHADE sst_clim[L=1] ! shade Quarter 1

This now provides an L=4 data set of seasonal climatologies centered in time on [15-Feb, 16-May, 16-Aug, 15-Nov].

So, while 10 years ago half the fun may have trying to be clever and code a faster solution than the last one that I had made, in my old age I tend to just want to get this portion of it done so that I can interpret the results. This is very nice indeed.

Technorati Tags:
, , ,