Home > paleoclimate > Huybers 2006b: Porting the Pleistocene

Huybers 2006b: Porting the Pleistocene

2011 January 25

NOAA’s Paleoclimate division maintains web and ftp sites for storing supplemental data and methods used in their research. You can find Huybers 2006b stuff here:
web: http://www.ncdc.noaa.gov/paleo/pubs/huybers2006b/huybers2006b.html
ftp: ftp://ftp.ncdc.noaa.gov/pub/data/paleo/climate_forcing/orbital_variations/huybers2006insolation/

Amid the data are two Matlab scripts and a Matlab data file.
daily_insolation.m
make_integrated_insolation_tables.m
orbital_parameter_data.mat (data)

This is the first time I’ve attempted to port a matlab script to R. It wasn’t that difficult. The “R.matlab” package allowed me to convert the matlab data set to a plain R text data file. The “signal” package gave me an “unwrap” routine. And an R “cheat sheet” for Matlab/Octave users included a little function to emulate the Matlab function meshgrid.

In addition, there was some reshaping of the function calls. Some restructuring of the conditionals. But pretty straightforward stuff when it comes to porting. I was pleased. Soon I had a set of data files in roughly the same format as the original. It was time to compare my data sets with Huybers’.

Oops.

And it wasn’t just the first set. It was all of them. Exactly the right shape. But with more variance in my set than the original.

A couple of possibilities occurred to me as I tried to debug the issue:

1. The R.matlab data conversion. I was able to rule this out when I found Berger’s 1991 data in the same ftp directory tree.

2a. The unwrap function. I had commented it out when I ran my numbers on the web server since I hadn’t loaded the R.matlab package into the web server. I’ll go into the web server set up more in another post. But, when I did do the package install and ran with the unwrap function enabled, there was little to no change in the data.

2b. The unwrap function may behave different in R than in Matlab. Still an open possibility.

3. 32bit -v- 64bit architecture. Probably not an issue. I compared R runs on a 32 and 64bit vm and saw no difference.

4. R precision -v- Matlab precision. Still an open issue.

So if any of you have porting experience and have seen something like this before, drop me a note.

Here are the R versions of Huybers’ code:
daily_insolation.R
make_integrated_insolation_tables.R
orbital_parameter_data.R (data)

Meanwhile, I just downloaded Octave (a more-or-less Matlab compatible software package) and am giving it a run. I’ll update this thread if I make any progress on the porting problem.

Advertisements