Re: [R-sig-Geo] Maptools package
On Thu, 4 Apr 2013, Jesse Berman wrote: Hi all, I've loaded the new version of R 3.0.0 and it seems the package 'maptools' is not available on the mirrors. It's not urgent, but I'm curious if anyone knows when it might become available? Thanks. This only applies to OSX for 3.0.0, 2.15.* appear still to be supported. The issue seems to be the absence of rgeos, which is suggested by maptools. There are also issues for OSX wrt. rgdal. This suggests that OSX users should stay on 2.15.* until things calm down. I released new sp, rgeos, and rgdal to suit an OSX transition from gcc to clang a few days ago, so probably OSX will catch up in a number of days. Hope this clarifies, Roger Regards, Jesse -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Maptools-package-tp7583198.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo -- Roger Bivand Department of Economics, NHH Norwegian School of Economics, Helleveien 30, N-5045 Bergen, Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43 e-mail: roger.biv...@nhh.no ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] seemingly unresolved problem with predict() in package raster
Hi all, I have a problem with the function raster::predict very similar to the one described here [1], [2] using raster package version 2.0-41 and party package version 1.0-6, where my model is a conditional inference forest (party::cforest). Could not find a solution in either post. The problem is the following: when I use the function predict() I get an error relating to a mismatch between the levels of factors in the data frame used to fit the model and those of the raster layers used to predict to the whole surface. The only difference that I can see is that the levels of the factors in the model are as follows (e.g.): morphrec.levels-levels(v$morphrec) morphrec.levels [1] 1 2 3 4 5 6 Whereas if I ask for the unique values of the same raster layer (where values were extracted from in an earlier step) that is now being used to predict a value for the response variable, I see the following: morphrec.values-sort(unique(getValues(morphrec))) morphrec.values [1] 1 2 3 4 5 6 When I try running predict, the following happens: r - predict(predictors, confor.dens, type='response',OOB=TRUE) Loading required package: tcltk Loading Tcl/Tk interface ... done Error in checkData(oldData, RET) : Classes of new data do not match original data predictors class : RasterStack dimensions : 2787, 2293, 6390591, 7 (nrow, ncol, ncell, nlayers) resolution : 200, 200 (x, y) extent : 296387.5, 754987.5, 7488413, 8045813 (xmin, xmax, ymin, ymax) coord. ref. : +proj=utm +zone=33 +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 names : bathy, slope, tri, grainsize, landscape, sedenv,morphrec min values : 0.44850001, 0.02432192, 0.04429007, 1., 21., 1., 1. max values : 2912.94849,46.13418,43.59974, 300.0, 431.0, 5.0, 6.0 confor.dens Random Forest using Conditional Inference Trees Number of trees: 1000 Response: density Inputs: bathy, slope, tri, grainsize, landscape, sedenv, morphrec Number of observations: 1020 Is there any way of telling predict() that the values of some raster layers will have to be considered characters? Or am I facing a different problem entirely? Any help will be much appreciated. I have been working on creating a small, reproducible example, but I can't seem to get the exact same behavior. Regards, Geno [1] http://r-sig-geo.2731867.n2.nabble.com/problems-with-predict-in-package-raster-td6627291.html#a6628041 [2] http://r-sig-geo.2731867.n2.nabble.com/new-function-predict-package-raster-td5796744.html#a5816769 Genoveva Gonzalez Mirelis, Scientist Institute of Marine Research Nordnesgaten 50 5005 Bergen, Norway +47 5523 6376 [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] predict function of regressors in the raster package
Hi, I would like to make a raster, based on the regression coefficients with 5 other rasters (Predictors). In theory this should be easy using the predict function of the raster package (http://cran.r-project.org/web/packages/raster/raster.pdf). But I fail to get it to work. First I fit a glm logit model using 5 cregressors: model - glm(c3p$t0~c3p$ycoord+c3p$prec+c3p$tempmean+c3p$tempmax+c3p$tempmin, family=binomial(link=logit)) logregt0 Call: glm(formula = c3p$t0 ~ c3p$ycoord + c3p$prec + c3p$tempmean + c3p$tempmax + c3p$tempmin, family = binomial(link = logit)) Coefficients: (Intercept)c3p$ycoord c3p$prec c3p$tempmean c3p$tempmax c3p$tempmin -7.179e+01 7.380e-06-1.621e-03-2.261e-01 7.630e-02 1.787e-01 Then I make a raster stack of the 5 prediction surfaces: newgrid = raster(newgrid.txt) Predictors - stack(newgrid) Predictors - addLayer(ycoord,prec,tempmean,tempmax,tempmin) layerNames(Predictors) [1] ycoord prec tempmean tempmax tempmin Then I try to make a new raster based on the reg. coeeficients of the 5 prediction surfaces: PredRegSurface - predict(object=Predictors, model=logregt0) Error in v[cells, ] - predv : number of items to replace is not a multiple of replacement length In addition: Warning message: 'newdata' had 4 rows but variable(s) found have 1350 rows In above error message: 1st part I do not understand, 2nd part refers to the Predictors raster (that has more than 40.000 rows) and the datafile on which the glm model is based v(which has 1350 rows). I must be doing something wrong in the syntax, which I have been changing for hours already without any result. Does anybody have any ideas or better even some script of this probably quite common thing (prediction from multiple regressors)? Hein Bouwmeester Wageningen University [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] collect latitude/longitude of various countries
Dear list, How can I obtain the latitude/longitude of various countries? In addition, does anyone know how to calculate absolute latitude? wanhai Best regards [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] rgdal and MSSQL Server geometries
Hi Shannon,
The following syntax has worked for me for the past year.
projstring - CRS('+init=epsg:28355')
##Establish the dsn
# Note: use the odbc tool in Windows to create a dsn for your SQL Server
database beforehand
myMSSQLdsn -
c(MSSQL:server=mysqlservername;database=mydatabase;trusted_connection=yes)
#Confirm connection is working
ogrListLayers(myMSSQLdsn)
#Reading sp object classes from SQL SERVER
lyr - c(mssqlserverTablename)
spdf -readOGR(dsn=myMSSQLdsn, layer=lyr,p4s=CRSargs(projstring))
#Writing sp object classes to SQL SERVER
lyrout= c(NewMSSQLSpatialTableName)
writeOGR(spdfname, dsn=myMSSQLdsn, layer=lyrout, driver=MSSQLSpatial,
layer_options=c(SRID=28355))
If you haven't created a dsn, this could be a big part of your problem.
Craig
--
View this message in context:
http://r-sig-geo.2731867.n2.nabble.com/rgdal-and-MSSQL-Server-geometries-tp7583193p7583203.html
Sent from the R-sig-geo mailing list archive at Nabble.com.
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] Help converting matlab raster to R raster brick
Hi All, I have been provided with a .mat file containing a time series of Sea Surface Temperature data (50 x 42 cells with 92 time layers). It was a Struct object in Matlab. I can happily import the file in to R using R.matlab, creating a list as follows. My question is how then to convert this to a raster brick? Thanks Craig library(R.matlab) test - readMat(CraigMundy.mat) str(test) List of 2 $ data :List of 9 ..$ : num [1, 1:92] 734138 734139 734140 734141 734142 ... ..$ : num [1:50, 1:42, 1:92] 14.4 14.5 15.2 15.3 14.8 ... ..$ : num [1:50, 1] -44.3 -44.3 -44.3 -44.2 -44.2 ... ..$ : num [1, 1:42] 147 147 147 147 147 ... ..$ : num [1, 1] 101 ..$ : chr [1, 1] SST: 3 Day Composite ..$ : num [1, 1] 0 ..$ : num [1, 1] 34.5 ..$ : num [1, 1] 0.44 ..- attr(*, dim)= int [1:3] 9 1 1 ..- attr(*, dimnames)=List of 3 .. ..$ : chr [1:9] Time Values Lat Long ... .. ..$ : NULL .. ..$ : NULL $ times: chr [1:92, 1] 31-Dec-2009 01-Jan-2010 02-Jan-2010 03-Jan-2010 ... - attr(*, header)=List of 3 ..$ description: chr MATLAB 5.0 MAT-file, Platform: PCWIN64, Created on: Thu Apr 04 11:29:53 2013 ..$ version: chr 5 ..$ endian : chr little version _ platform x86_64-w64-mingw32 arch x86_64 os mingw32 system x86_64, mingw32 status major 3 minor 0.0 year 2013 month 04 day03 svn rev62481 language R version.string R version 3.0.0 (2013-04-03) nickname Masked Marvel -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Help-converting-matlab-raster-to-R-raster-brick-tp7583204.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] predict function of regressors in the raster package
Hi Hein, I'm not sure if this will help, but one thing to check is that your prediction grid has covariate data for each of the 40,000 cells. If a large number of cells have 'NA' as data values, then sometimes the prediction will not work. Offhandedly, it strikes me that ycoord may be limited to only your t0 locations. Jesse -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/predict-function-of-regressors-in-the-raster-package-tp7583201p7583205.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] predict function of regressors in the raster package
Hein, I think this problem goes away if you use a clean formula, like this: model - glm(t0~ycoord+prec+tempmean+tempmax+tempmin, family=binomial(link=logit), data=c3p) Robert On Thu, Apr 4, 2013 at 7:42 AM, Jesse Berman berman.je...@gmail.com wrote: Hi Hein, I'm not sure if this will help, but one thing to check is that your prediction grid has covariate data for each of the 40,000 cells. If a large number of cells have 'NA' as data values, then sometimes the prediction will not work. Offhandedly, it strikes me that ycoord may be limited to only your t0 locations. Jesse -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/predict-function-of-regressors-in-the-raster-package-tp7583201p7583205.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] seemingly unresolved problem with predict() in package raster
cforest with factors are currently not supported (although they may work in some cases). I will change the predict function to make it more general by adding a 'levels' argument such that you can indicate (for models with a non-standard structure such as cforest) which variables are factors and what their levels are. Robert On Thu, Apr 4, 2013 at 4:22 AM, Gonzalez-Mirelis Genoveva genoveva.gonzalez-mire...@imr.no wrote: Hi all, I have a problem with the function raster::predict very similar to the one described here [1], [2] using raster package version 2.0-41 and party package version 1.0-6, where my model is a conditional inference forest (party::cforest). Could not find a solution in either post. The problem is the following: when I use the function predict() I get an error relating to a mismatch between the levels of factors in the data frame used to fit the model and those of the raster layers used to predict to the whole surface. The only difference that I can see is that the levels of the factors in the model are as follows (e.g.): morphrec.levels-levels(v$morphrec) morphrec.levels [1] 1 2 3 4 5 6 Whereas if I ask for the unique values of the same raster layer (where values were extracted from in an earlier step) that is now being used to predict a value for the response variable, I see the following: morphrec.values-sort(unique(getValues(morphrec))) morphrec.values [1] 1 2 3 4 5 6 When I try running predict, the following happens: r - predict(predictors, confor.dens, type='response',OOB=TRUE) Loading required package: tcltk Loading Tcl/Tk interface ... done Error in checkData(oldData, RET) : Classes of new data do not match original data predictors class : RasterStack dimensions : 2787, 2293, 6390591, 7 (nrow, ncol, ncell, nlayers) resolution : 200, 200 (x, y) extent : 296387.5, 754987.5, 7488413, 8045813 (xmin, xmax, ymin, ymax) coord. ref. : +proj=utm +zone=33 +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 names : bathy, slope, tri, grainsize, landscape, sedenv,morphrec min values : 0.44850001, 0.02432192, 0.04429007, 1., 21., 1., 1. max values : 2912.94849,46.13418,43.59974, 300.0, 431.0, 5.0, 6.0 confor.dens Random Forest using Conditional Inference Trees Number of trees: 1000 Response: density Inputs: bathy, slope, tri, grainsize, landscape, sedenv, morphrec Number of observations: 1020 Is there any way of telling predict() that the values of some raster layers will have to be considered characters? Or am I facing a different problem entirely? Any help will be much appreciated. I have been working on creating a small, reproducible example, but I can't seem to get the exact same behavior. Regards, Geno [1] http://r-sig-geo.2731867.n2.nabble.com/problems-with-predict-in-package-raster-td6627291.html#a6628041 [2] http://r-sig-geo.2731867.n2.nabble.com/new-function-predict-package-raster-td5796744.html#a5816769 Genoveva Gonzalez Mirelis, Scientist Institute of Marine Research Nordnesgaten 50 5005 Bergen, Norway +47 5523 6376 [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] rgdal and MSSQL Server geometries
Craig, thanks for the idea. I did create a DSN using ODBC as well as using an
explicit string connection. Both give me the same result of: 'Cannot open data
source'
Thus, I am thoroughly confused because I can connect to the SQL Server via all
of my other methods (RODBC) and software with no issues. I am wondering if
there is something with rgdal that is having a miscommunication with the proper
drivers. If so, I really do not know how to solve the issue. Hopefully some
other ideas will surface from the r-sig-geo list that can provide additional
troubleshooting that I can try.
Cheers,
Shannon
-Original Message-
From: r-sig-geo-boun...@r-project.org [mailto:r-sig-geo-boun...@r-project.org]
On Behalf Of cmundy
Sent: Thursday, April 04, 2013 7:52 AM
To: r-sig-geo@r-project.org
Subject: Re: [R-sig-Geo] rgdal and MSSQL Server geometries
Hi Shannon,
The following syntax has worked for me for the past year.
projstring - CRS('+init=epsg:28355')
##Establish the dsn
# Note: use the odbc tool in Windows to create a dsn for your SQL Server
database beforehand myMSSQLdsn -
c(MSSQL:server=mysqlservername;database=mydatabase;trusted_connection=yes)
#Confirm connection is working
ogrListLayers(myMSSQLdsn)
#Reading sp object classes from SQL SERVER lyr - c(mssqlserverTablename)
spdf -readOGR(dsn=myMSSQLdsn, layer=lyr,p4s=CRSargs(projstring))
#Writing sp object classes to SQL SERVER lyrout= c(NewMSSQLSpatialTableName)
writeOGR(spdfname, dsn=myMSSQLdsn, layer=lyrout, driver=MSSQLSpatial,
layer_options=c(SRID=28355))
If you haven't created a dsn, this could be a big part of your problem.
Craig
--
View this message in context:
http://r-sig-geo.2731867.n2.nabble.com/rgdal-and-MSSQL-Server-geometries-tp7583193p7583203.html
Sent from the R-sig-geo mailing list archive at Nabble.com.
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] rgdal and MSSQL Server geometries
Hi Shannon, The other requirement to read tables from MS SQL Server is that details of the table with geometry columns must be registered in a table called geometry_columns in your database. If your database doesn't have this table, do the following; 1) write a small point or polygon object back to SQL Server using writeOGR. 2) This process will create two additional tables, the geometry_columns table and a spatial_ref_sys table. 3) open the geometry_columns table in SQL Server, and create a new record with the details of the table you want to read into R, following the pattern of the record created by the write process. 4) try reading the table again I'm out of the office for a day or so, but happy to help work this out. I use rgdal to exchange spatial data between R and SQL Server on a weekly basis, so can at least confirm that it works very nicely. Craig -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/rgdal-and-MSSQL-Server-geometries-tp7583193p7583209.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Re: [R-sig-Geo] Help converting matlab raster to R raster brick
Hi Craig, see the array method for ?brick: brick(x, xmn=0, xmx=1, ymn=0, ymx=1, crs=NA So, something like x - brick(test$data[[2], xmn = min(test$data[[5]]), ymx = max(test$data[[5]]), ymn = min(test$data[[4]]), ymx = max(test$data[[4]])) You might need to mess around to get the x/y ranges right, and figure out if you need to transpose or do other reorientations. This assumes the lon/lat vectors are regular so you should check that. Cheers, Mike On Thu, Apr 4, 2013 at 9:00 AM, cmundy craig.mu...@utas.edu.au wrote: Hi All, I have been provided with a .mat file containing a time series of Sea Surface Temperature data (50 x 42 cells with 92 time layers). It was a Struct object in Matlab. I can happily import the file in to R using R.matlab, creating a list as follows. My question is how then to convert this to a raster brick? Thanks Craig library(R.matlab) test - readMat(CraigMundy.mat) str(test) List of 2 $ data :List of 9 ..$ : num [1, 1:92] 734138 734139 734140 734141 734142 ... ..$ : num [1:50, 1:42, 1:92] 14.4 14.5 15.2 15.3 14.8 ... ..$ : num [1:50, 1] -44.3 -44.3 -44.3 -44.2 -44.2 ... ..$ : num [1, 1:42] 147 147 147 147 147 ... ..$ : num [1, 1] 101 ..$ : chr [1, 1] SST: 3 Day Composite ..$ : num [1, 1] 0 ..$ : num [1, 1] 34.5 ..$ : num [1, 1] 0.44 ..- attr(*, dim)= int [1:3] 9 1 1 ..- attr(*, dimnames)=List of 3 .. ..$ : chr [1:9] Time Values Lat Long ... .. ..$ : NULL .. ..$ : NULL $ times: chr [1:92, 1] 31-Dec-2009 01-Jan-2010 02-Jan-2010 03-Jan-2010 ... - attr(*, header)=List of 3 ..$ description: chr MATLAB 5.0 MAT-file, Platform: PCWIN64, Created on: Thu Apr 04 11:29:53 2013 ..$ version: chr 5 ..$ endian : chr little version _ platform x86_64-w64-mingw32 arch x86_64 os mingw32 system x86_64, mingw32 status major 3 minor 0.0 year 2013 month 04 day03 svn rev62481 language R version.string R version 3.0.0 (2013-04-03) nickname Masked Marvel -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Help-converting-matlab-raster-to-R-raster-brick-tp7583204.html Sent from the R-sig-geo mailing list archive at Nabble.com. ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo -- Michael Sumner Hobart, Australia e-mail: mdsum...@gmail.com [[alternative HTML version deleted]] ___ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
[R-sig-Geo] [raster] conserving dimensions when regridding a 4D Brick or Stack
summary: When regridding a 4D Brick (e.g., the netCDF data variable
Fraction_of_Emissions(TSTEP, LAY, ROW, COL) in
netcdf GFED-3.1_2008_N2O_3hourly_fractions {
dimensions:
TSTEP = 12 ;
LAY = 8 ;
ROW = 360 ;
COL = 720 ;
variables:
float Fraction_of_Emissions(TSTEP, LAY, ROW, COL) ;
Fraction_of_Emissions:units = unitless ;
Fraction_of_Emissions:long_name = 3hourly fraction of monthly N2O
emission ;
Fraction_of_Emissions:_FillValue = 9.96921e+36f ;
int TSTEP(TSTEP) ;
TSTEP:long_name = index of month in 2008 ;
TSTEP:units = month ;
int LAY(LAY) ;
LAY:units = 4-digit hour ;
LAY:long_name = time of day starting 3-hour-long period ;
float ROW(ROW) ;
ROW:units = degrees_north ;
ROW:long_name = latitude ;
float COL(COL) ;
COL:units = degrees_east ;
COL:long_name = longitude ;
) I want to keep all 4 dimensions in the output. But I am only able to
write 3D output, since I lose dimension=LAY:
netcdf GFED-3.1_2008_N2O_3hourly_fractions_regrid {
dimensions:
COL = 459 ;
ROW = 299 ;
TSTEP = UNLIMITED ; // (12 currently)
variables:
double COL(COL) ;
COL:units = meter ;
COL:long_name = COL ;
double ROW(ROW) ;
ROW:units = meter ;
ROW:long_name = ROW ;
int TSTEP(TSTEP) ;
TSTEP:units = month ;
TSTEP:long_name = TSTEP ;
float Fraction_of_Emissions(TSTEP, ROW, COL) ;
Fraction_of_Emissions:units = unitless ;
Fraction_of_Emissions:_FillValue = -3.4e+38 ;
Fraction_of_Emissions:missing_value = -3.4e+38 ;
Fraction_of_Emissions:long_name = 3hourly fraction of monthly N2O
emission ;
Fraction_of_Emissions:projection = +proj=lcc +lat_1=33 +lat_2=45
+lat_0=40 +lon_0=-97 +a=637 +b=637 ;
Fraction_of_Emissions:projection_format = PROJ.4 ;
) Can one
+ create a 4D RasterBrick or RasterStack
+ instruct projectRaster
so as to get 4D regridded output? Or must one
- create 3D RasterLayer's (i.e., one per LAY) from the 4D Brick
- regrid each 3D Layer
- reassemble the 4D Brick
? or Something Completely Different?
details:
As part of
https://bitbucket.org/tlroche/gfed-3.1_global_to_aqmeii-na
I need to regrid all 4 dimensions of the data variable
Fraction_of_Emissions(TSTEP, LAY, ROW, COL) in a netCDF file:
netcdf GFED-3.1_2008_N2O_3hourly_fractions {
dimensions:
TSTEP = 12 ;
LAY = 8 ;
ROW = 360 ;
COL = 720 ;
variables:
float Fraction_of_Emissions(TSTEP, LAY, ROW, COL) ;
Fraction_of_Emissions:units = unitless ;
Fraction_of_Emissions:long_name = 3hourly fraction of monthly N2O
emission ;
Fraction_of_Emissions:_FillValue = 9.96921e+36f ;
int TSTEP(TSTEP) ;
TSTEP:long_name = index of month in 2008 ;
TSTEP:units = month ;
int LAY(LAY) ;
LAY:units = 4-digit hour ;
LAY:long_name = time of day starting 3-hour-long period ;
float ROW(ROW) ;
ROW:units = degrees_north ;
ROW:long_name = latitude ;
float COL(COL) ;
COL:units = degrees_east ;
COL:long_name = longitude ;
In
https://bitbucket.org/tlroche/gfed-3.1_global_to_aqmeii-na/src/6a214b2076dee1ce4cae5ee44779642b33481d81/vis_regrid_vis.r?at=master
I load the data like (substituting and reformatting for mail)
hour3ly.in.raster - raster::brick(
'GFED-3.1_2008_N2O_3hourly_fractions.nc',
varname='Fraction_of_Emissions')
hour3ly.in.raster@crs - global.crs
hour3ly.in.raster
# class : RasterBrick
# dimensions : 360, 720, 259200, 12 (nrow, ncol, ncell, nlayers)
# resolution : 0.5, 0.5 (x, y)
# extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=longlat +ellps=WGS84
# data source : .../GFED-3.1_2008_N2O_3hourly_fractions.nc
# names : X1, X2, X3, X4, X5, X6, X7, X8, X9, X10, X11, X12
# month : 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12
# varname : Fraction_of_Emissions
# level : 1
and then regrid with
hour3ly.out.raster -
raster::projectRaster(
from=hour3ly.in.raster, to=template.raster, crs=out.crs,
method='bilinear', overwrite=TRUE, format='CDF',
# args from writeRaster
varname=hour3ly_out_datavar_name,
varunit=hour3ly_out_datavar_units,
longname=hour3ly_out_datavar_long_name,
xname=hour3ly_out_x_var_name,
yname=hour3ly_out_y_var_name,
zname=hour3ly_out_z_var_name,
zunit=hour3ly_out_z_var_units,
filename=hour3ly_out_fp)
hour3ly.out.raster
# class : RasterBrick
# dimensions : 299, 459, 137241, 12 (nrow, ncol, ncell, nlayers)
# resolution : 12000, 12000 (x, y)
# extent : -2556000, 2952000, -1728000, 186 (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=lcc +lat_1=33 +lat_2=45 +lat_0=40 +lon_0=-97 +a=637
+b=637
# data source :
/tmp/gfed-3.1_global_to_aqmeii-na/GFED-3.1_2008_N2O_3hourly_fractions_regrid.nc
# names : X1, X2, X3, X4, X5, X6, X7, X8, X9, X10, X11, X12
# month : 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12
# varname : Fraction_of_Emissions
The output appears to
