R help - Jan 2013 - Netcdf and Raster Package Questions, Need .asc File for GIS

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
---------------------------------

Hello,

I don't understand why you need to use R first, as GRASS and ArcGIS are 
able to read NetCDF files.

Regards,
Pascal


Le 30/01/2013 00:37, Douglas M. Hultstrand a ?crit :
> 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
>
>
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

Douglas M. Hultstrand <dmhultst <at> metstat.com> writes:
> 
> 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.
Please post on R-sig-geo, where the raster developers often handle questions
like this.

Roger