[R] Netcdf and Raster Package Questions, Need .asc File for GIS

[Douglas M. Hultstrand]

Hello R-Group,

I am new working with netcdf files and the raster package in R.I am 
trying to read in a netcdf file using the package "ncdf".I am able to 
get the lat, lon and parameter I need and can plot using 
fill.contour.Ultimately, I am trying to create a .asc file to reafd into 
GIS.I am using the package "raster" to read the parameter.When I read in 
with "raster", the extent of the raster is between 0 and 1 for both x 
and y directions.Also, I have to transpose the grid so it is oriented 
the correctly.In order to create the .asc file I have to resample to 
square grids for "writeRaster" to work.

Below I supplied my objective, questions and R code used.(I also 
attached .docx file with code and images for reference)

*_Objective_***- read netcdf file and create .asc file to read into GIS

*_Questions:_*

1) using the raster package how can I set the projection, extent, and 
cell size properly.Here are the dimensions and global attributes of the 
netcdf file, I used "ncdump -h "

dimensions:

lon = 4200 ;

lat = 2400 ;

time = UNLIMITED ; // (1 currently)

rainf_profiles = 1 ;

tair_profiles = 1 ;

global attributes:

:missing_value = -9999.f ;

:TITLE = "LIS land surface model output" ;

:MAP_PROJECTION = "EQUIDISTANT CYLINDRICAL" ;

:SOUTH_WEST_CORNER_LAT = 48.005f ;

:SOUTH_WEST_CORNER_LON = -167.995f ;

:DX = 0.01f ;

:DY = 0.01f ;

2) Is there a better way to read a netcdf file and create an .asc file 
that can be read into GIS (GRASS or ArcGIS)



*_Code Used_*

# get net cdf file, big file

wget 
ftp://ftp.nohrsc.nws.gov/pub/staff/fall/Hultstrand/NOAH32.201204110000.d01.nc

##################

#### COMMANDS ##

##################

# look at netcdf file header info

ncdump -h NOAH32.201204110000.d01.nc

##########

# Start R #

##########

library(ncdf)

ex.nc = open.ncdf("NOAH32.201204110000.d01.nc")

print(ex.nc)

summary(ex.nc)

x = get.var.ncdf(ex.nc, "lon")

y = get.var.ncdf(ex.nc, "lat")

swe = get.var.ncdf(ex.nc, "SWE") # kg/m2

swe_in <- swe * 0.0393701

# image plot of swe, takes time to load

filled.contour(x,y,swe_in, color = terrain.colors, asp = 1)

# load raster library

library(raster)

r = raster(swe_in)

# raster is sideways, rotate

m <- t(swe_in)

m <- m[nrow(m):1,]

r <- raster(m)

plot(r)

writeRaster(r,file="test_grid2.asc", format="ascii") # need square grids

# resample to square grids

r2 = raster(r)

res(r2) = min(res(r))

res(r2)

r2 <- resample(r, r2, method='bilinear')

writeRaster(r2,file="test_grid24.asc", format="ascii")


Thanks,

Doug

-- 
---------------------------------
Douglas M. Hultstrand, MS
Senior Hydrometeorologist
Metstat, Inc. Windsor, Colorado
voice: 720.771.5840
email: dmhultst at metstat.com
web: http://www.metstat.com
---------------------------------