create ascii grid from grid.list object (fields)?

Dear list, I've used R for simple statistics in the past, but this is my 
first attempt at manipulating spatial data. The transition is 
frustrating and I'd appreciate any advice you can give me on the problem 
below.

I'm creating interpolated surfaces of PCA loadings extracted from some 
precipitation data in Oregon/WA USA. I'm using the package fields. I can 
create the interpolated surface, but I need to use the data from the 
interpolated surface, either extracting values from the surface for 
specific x,y locations, or exporting the surface to an ascii grid file, 
or some format that I can either access via another R package (rgdal??) 
or ArcGIS, which I've used in the past.

I'd be very grateful for suggestions and/or a snippet of example code 
telling me how to make the interpolated surface created below, into a 
spatial object that I can analyze further.

Thanks in advance,

Zack

Another time, please use a cluefull text-only email client - the code was 
completely mixed up.

#######################################################
require(fields)
require(maptools)
require(sp)

setwd("C:/climate_data/pnw_data")
data <- readShapePoints("climate_PNW_PCLOAD.shp")

If possible use the rgdal package, it gives many more drivers for vector 
and raster data - see readOGR().

# read in shapefile with X,Y and PC loadings for climate 
stationsproj4string=CRS(as.character(NA)), verbose=FALSE)

No idea what this was supposed to be

# data <- as.data.frame(data)
# long <- x[,1]
# longitude column
# lat <- x[,2]
# latitude column
# minlat <- min(lat)
# maxlat <- max(lat)
# minlong <- min(long)
# maxlong <- max(long)
# grid.l<- list( X1= seq(minlong, maxlong,,100), 
# X3=seq(minlat,maxlat,,100))

PredGrid <- Sobj_SpatialGrid(data)$SG
grd <- slot(PredGrid, "grid")
o <- rep(NA, prod(slot(grd, "cells.dim")))
PredGridDF <- SpatialGridDataFrame(slot(PredGrid, "grid"),
  data=data.frame(o=o))

This makes an empty SpatialGridDataFrame, use the maxDim= argument to 
Sobj_SpatialGrid() to control resolution.

# create grid.list object using ranges of Long. and Lat.
# make.surface.grid(grid.l)
# y <- data$pc1
# Princ. Comp. loading I want to interpolate
# LL <- cbind(long,lat)
# combine Lat and Long
# fit = Krig(LL,y, theta = 10)

fit <- Krig(coordinates(data), data$pc1, theta=10)

Use the data object directly

# fit Princ. comp. 1 to Lat long.
# surface( fit, grid.list=grid.l)

PredGridDF$pc1 <- predict(fit, coordinates(PredGridDF))
spplot(PredGridDF, "pc1", col.regions=tim.colors())

Use the predict() method on the prediction grid coordinates, assigning to 
a variable in the object (treat it as you would treat a data frame).

writeGDAL(PredGridDF["pc1"], "pc1.tif", drivername="GTiff")

to write a geotiff - see other drivers with gdalDrivers().

Hope this helps,

Roger

#  create interpolated surface based on Krig object

####################################################
	[[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo at stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no