Code will be in the drop box in a bit, once I shower:  [DONE] This is a wholesale replacement of previous versions, completely rewritten in raster. It will be the base going forward. All of the analysis routines will be rewritten using raster.  For time series functionality I will continue to use zoo as that package rocks and Gabor is a great help. Also, Thanks to the following: Robert H for a great effort on raster, Ron Broberg for his careful eye and GIS inspiration, Zeke for verifying my early versions so long ago, SteveMc for encouragement to learn R and for code here and there, RomanM for his patience and R help, JeffId for pieces of code here and there and some nice R code that helped me understand, Ryan O for nice bits of code, Peter Oneill for KML code, ( needs to get re integrated), Nick Stokes for a good deal of R code that  helped me early on, and too many people to mention on the R help list> If I left you off, ping me and I’ll edit you in. oh and Nicholas for fixing the bug in “uncompress”

Download the Zip: unzip. Make the folder your working directory. Run Downloadall.R, Run SetUpData.R

Run LandOcean.

Report any issues.

When I put Moshtemp to bed last night there was a lingering problem. I was significantly warmer than CRU. About .2C warmer. Since we have different data as Ron Broberg pointed out I took some comfort in that. But our methods are nearly identical. Our SST matched. My Land while noisier than theirs matched when the mean of monthly values was taken. The noise in my data cancelled out. So why was my global figure higher.  I spent the day cleaning up code and improving the style, basically getting a consistent coding style with respect to naming. I lean toward this guide and suggest it for others

https://docs.google.com/View?id=dddzqd53_2646dcw759cb

henrik write very readable code and I use his packages, so I’m moving toward his convention. I don’t follow it entirely, but this version is much improved. The amount of code I was able to remove by moving to raster was substantial and many of my utility programs went into the trash can. I’ve added some new general purpose code for IO that is part of the ICOADS project ( what a bear )

Let me dwell a minute on the bug. Here is the wrong code:

Land<-Weights*LandGrid*landMask

Weights <- area(SST,na.rm=T,weight=T)

Ocean <- SST*Weights*OceanPercentMask

Coastal <-Ocean+Land

Global <- cover(Coastal, Ocean, Land)

The problem with that is in the calculation of weights. The weights cannot be computed separately and then added. Instead we must proceed like this:

Land <- pointsToRaster(landMask,xy=getInvLonLat(Inv),values=Anom,fun=mean,na.rm=T)

Coastal <- (Land*landMask) + (SST*OceanPercentMask)

Global  <- cover(Coastal,Land*landMask,SST*OceanPercentMask)

Global  <- Global * area(Global,na.rm=T,weight=T)

Land    <- Land   * area(Land,na.rm=T,weight=T)

SST     <-   SST * area(SST,na.rm=T,weight=T)

Very simply. To create a Coastal cell I want to add the Anomaly in the land cell with the anomaly in the SST cell. But I want to weight them ONLY by the land percentage/ocean percentage

When I create a Global brick I want to “cover” this coastal brick, with an Land brick that is weighted by percent of land and by a SST brick that is weighted by percent of ocean of ocean in each cell.

Then and only then do I want to weight each cell by the area. Why? because we care about the total area sampled and we want these weights to equal one. In the “bugged” version the weights for SST and Land were computed seperately and thats wrong.

Finally I take the land brick and the SST brick and I weight them separately for the individual comparisons. This has to be done after a global brick is constructed. the land weights have to equal 1 and the SST weights have to equal one.

Time for charts:

Land with Moshtemp  in BLUE

As we’ve noted before we run a bit warmer than CRU.  next is the SST

And Finally Global:

And  a sample of global maps ( yes the colors and scales need work )