Hi Frede,

Thanks so much! It seems to work perfectly, exactly what I wanted to do! :)

All the best,

Bea

On 16.04.2014 12:25, Frede Aakmann Tøgersen wrote:
Hi Bea

Well the first hit lead me to rasterProject{raster}. Will this suit you?


rasterMG.proj <- projectRaster(rasterMG, crs=CRS("+init=epsg:21781"))
print(rasterMG.proj)
class       : RasterLayer
dimensions  : 116, 91, 10556  (nrow, ncol, ncell)
resolution  : 40.1, 40.1  (x, y)
extent      : 478794.9, 482444, 646645.8, 651297.4  (xmin, xmax, ymin, ymax)
coord. ref. : +init=epsg:21781 +proj=somerc +lat_0=46.95240555555556 
+lon_0=7.439583333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel 
+towgs84=674.374,15.056,405.346,0,0,0,0 +units=m +no_defs
data source : in memory
names       : part.a
values      : 0, 1  (min, max)


The difference can be seen by plotting.

plot(rasterMG)
plot(rasterMG.proj)




Yours sincerely / Med venlig hilsen


Frede Aakmann Tøgersen
Specialist, M.Sc., Ph.D.
Plant Performance & Modeling

Technology & Service Solutions
T +45 9730 5135
M +45 2547 6050
fr...@vestas.com
http://www.vestas.com

Company reg. name: Vestas Wind Systems A/S
This e-mail is subject to our e-mail disclaimer statement.
Please refer to www.vestas.com/legal/notice
If you have received this e-mail in error please contact the sender.


-----Original Message-----
From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
On Behalf Of Frede Aakmann Tøgersen
Sent: 16. april 2014 12:11
To: Bea GD; r-help@r-project.org
Subject: Re: [R] netCDF to raster and spatial projection

Well no need to have your data. Usually one can find some suitable data in
the help for the functions under question. So here is a reproducible example
from ?meuse.grid.

data(meuse.grid)
coordinates(meuse.grid) <- ~x+y
proj4string(meuse.grid) <- CRS("+init=epsg:28992") # see ?meuse for this
crs
gridded(meuse.grid) <- TRUE
summary(meuse.grid)
Object of class SpatialPixelsDataFrame
Coordinates:
      min    max
x 178440 181560
y 329600 333760
Is projected: TRUE
proj4string :
[+init=epsg:28992 +proj=sterea +lat_0=52.15616055555555
+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000
+ellps=bessel +units=m +no_defs]
Number of points: 3103
Grid attributes:
   cellcentre.offset cellsize cells.dim
x            178460       40        78
y            329620       40       104
Data attributes:
      part.a           part.b            dist        soil     ffreq
  Min.   :0.0000   Min.   :0.0000   Min.   :0.0000   1:1665   1: 779
  1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:0.1193   2:1084   2:1335
  Median :0.0000   Median :1.0000   Median :0.2715   3: 354   3: 989
  Mean   :0.3986   Mean   :0.6014   Mean   :0.2971
  3rd Qu.:1.0000   3rd Qu.:1.0000   3rd Qu.:0.4402
  Max.   :1.0000   Max.   :1.0000   Max.   :0.9926

rasterMG <- raster(meuse.grid)
print(rasterMG)
class       : RasterLayer
dimensions  : 104, 78, 8112  (nrow, ncol, ncell)
resolution  : 40, 40  (x, y)
extent      : 178440, 181560, 329600, 333760  (xmin, xmax, ymin, ymax)
coord. ref. : +init=epsg:28992 +proj=sterea +lat_0=52.15616055555555
+lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000
+ellps=bessel +units=m +no_defs
data source : in memory
names       : part.a
values      : 0, 1  (min, max)

rasterMG.proj <- spTransform(rasterMG, CRS("+init=epsg:21781"))
Error in spTransform(rasterMG, CRS("+init=epsg:21781")) :
   load package rgdal for spTransform methods
Now perhaps doing the projection before casting it to a raster will work (yes
I'm guessing since error thrown above is not very informative and traceback()
does not give anything useful).

meuse.grid.proj <- spTransform(meuse.grid, CRS("+init=epsg:21781"))
Warning message:
In spTransform(meuse.grid, CRS("+init=epsg:21781")) :
   Grid warping not available, coercing to points

summary(meuse.grid.proj)
Object of class SpatialPointsDataFrame
Coordinates:
        min      max
x 479029.8 482197.1
y 646927.4 651003.9
Is projected: TRUE
proj4string :
[+init=epsg:21781 +proj=somerc +lat_0=46.95240555555556
+lon_0=7.439583333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel
+towgs84=674.374,15.056,405.346,0,0,0,0 +units=m +no_defs]
Number of points: 3103
Data attributes:
      part.a           part.b            dist        soil     ffreq
  Min.   :0.0000   Min.   :0.0000   Min.   :0.0000   1:1665   1: 779
  1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:0.1193   2:1084   2:1335
  Median :0.0000   Median :1.0000   Median :0.2715   3: 354   3: 989
  Mean   :0.3986   Mean   :0.6014   Mean   :0.2971
  3rd Qu.:1.0000   3rd Qu.:1.0000   3rd Qu.:0.4402
  Max.   :1.0000   Max.   :1.0000   Max.   :0.9926

gridded(meuse.grid.proj) <- TRUE
suggested tolerance minimum: 0.737332
Error in points2grid(points, tolerance, round) :
   dimension 1 : coordinate intervals are not constant
A warning and an error indicates that the projection results in something that
is not on a regular grid.

I don't know what to do but to read some more of the documentation for sp,
rgdal, etc.

Hopefully somebody comes up with something that helps us.

Yours sincerely / Med venlig hilsen


Frede Aakmann Tøgersen
Specialist, M.Sc., Ph.D.
Plant Performance & Modeling

Technology & Service Solutions
T +45 9730 5135
M +45 2547 6050
fr...@vestas.com<mailto:fr...@vestas.com>
http://www.vestas.com<http://www.vestas.com/>

Company reg. name: Vestas Wind Systems A/S
This e-mail is subject to our e-mail disclaimer statement.
Please refer to
www.vestas.com/legal/notice<http://www.vestas.com/legal/notice>
If you have received this e-mail in error please contact the sender.

From: Bea GD [mailto:aguitatie...@hotmail.com]
Sent: 16. april 2014 11:00
To: Frede Aakmann Tøgersen; r-help@r-project.org
Subject: Re: [R] netCDF to raster and spatial projection

Hi Frede,

Thanks so much for your reply!

I've tried what you said but I get the following error:



Error in spTransform(rasterDF1, crs("+init=epsg:21781")) :

   load package rgdal for spTransform methods

I've checked search() and rgdal is before sp.


search()
  [1] ".GlobalEnv"        "package:chron"     "package:sm"        
"package:rgeos"

  [5] "package:maptools"  "package:ncdf"      "package:rgdal"
"package:raster"

  [9] "package:sp"        "tools:rstudio"     "package:stats"     
"package:graphics"

[13] "package:grDevices" "package:utils"     "package:datasets"
"package:methods"

[17] "Autoloads"         "package:base"

Also, when I do library("rgdal") I get this message:


library("rgdal")
rgdal: version: 0.8-16, (SVN revision 498)

Geospatial Data Abstraction Library extensions to R successfully loaded

Loaded GDAL runtime: GDAL 1.10.1, released 2013/08/26

Path to GDAL shared files: C:/Users/bgonzale/Documents/R/win-
library/3.0/rgdal/gdal

GDAL does not use iconv for recoding strings.

Loaded PROJ.4 runtime: Rel. 4.8.0, 6 March 2012, [PJ_VERSION: 480]

Path to PROJ.4 shared files: C:/Users/bgonzale/Documents/R/win-
library/3.0/rgdal/proj

Maybe the problem is something to do with this? I don't know how to fix it.

I'd like to provide you with some data, how would be the best way to
post/sahre a raster?

Thanks!


On 16.04.2014 10:17, Frede Aakmann Tøgersen wrote:

Hi Beatriz



Try to skip this step



      # Reprojecting into CH1903_LV03

      # First, change the coordinate reference system (crs)

      proj4string(rasterDF1) <- "+init=epsg:21781"





And just do this



      # Second, reproject the raster

      rasterDF1.CH <- spTransform(rasterDF1, crs("+init=epsg:21781"))





Also there is a spTransform both in the rgdal and sp packages. So are they
masking each other? rgdal should be before sp in the search() list.



I cannot be of more help since you provided no data.





Yours sincerely / Med venlig hilsen





Frede Aakmann Tøgersen

Specialist, M.Sc., Ph.D.

Plant Performance & Modeling



Technology & Service Solutions

T +45 9730 5135

M +45 2547 6050

fr...@vestas.com<mailto:fr...@vestas.com>

http://www.vestas.com



Company reg. name: Vestas Wind Systems A/S

This e-mail is subject to our e-mail disclaimer statement.

Please refer to
www.vestas.com/legal/notice<http://www.vestas.com/legal/notice>

If you have received this e-mail in error please contact the sender.



-----Original Message-----

From: r-help-boun...@r-project.org<mailto:r-help-boun...@r-project.org>
[mailto:r-help-boun...@r-project.org]

On Behalf Of Beatriz R. Gonzalez Dominguez

Sent: 16. april 2014 08:22

To: r-help@r-project.org<mailto:r-help@r-project.org>

Subject: [R] netCDF to raster and spatial projection



I've recently started using R for spatial data. I'd be really grateful

if you could could help me with this. Thanks!



Sorry I don't provide a reproducible example. Please ask me if you have

any questions about the data.



I've extracted data from a multidimensional netCDF file. This file had

longitude, latitude and temperature data (for 12 months of a specific year).



  From this netCDF I've got a data frame for January with these

variables: longitude, latitude, temperature.



With this data frame I've created a raster.



      # Packages

      library("sp")

      library("raster")

      library("rgdal")

      library("ncdf")

      library("maptools")

      library("rgeos")

      library("sm")

      library("chron")



      # Dataframe to raster

      # Create spatial points data frame

      coordinates(tmp.df01) <- ~ lon + lat

      # Coerce to SpatialPixelsDataFrame

      gridded(tmp.df01) <- T

      # Coerce to raster

      rasterDF1 <- raster(tmp.df01)

      > print(tmp.df01)

      class       : SpatialPixelsDataFrame

      dimensions  : 103, 241, 24823, 1  (nrow, ncol, npixels, nlayers)

      resolution  : 0.02083333, 0.02083333  (x, y)

      extent      : 5.739583, 10.76042, 45.73958, 47.88542  (xmin, xmax,

ymin, ymax)

      coord. ref. : NA

      names       :           TabsM_1

      min values  : -18.1389980316162

      max values  :  2.26920962333679



There is no value for 'coord. ref.'



The projection of the original netCDF was WGS84. So I gave this

projection to the raster.



      proj4string(rasterDF1) <- "+proj=longlat +datum=WGS84 +ellps=WGS84

+towgs84=0,0,0"



Then, I wanted to reproject my raster to another projection:



      # Reprojecting into CH1903_LV03

      # First, change the coordinate reference system (crs)

      proj4string(rasterDF1) <- "+init=epsg:21781"

      # Second, reproject the raster

      rasterDF1.CH <- spTransform(rasterDF1, crs("+init=epsg:21781"))



At this point I get the following error:



      Error in spTransform(rasterDF1, crs("+init=epsg:21781")) :

        load package rgdal for spTransform methods



But the package rgdal is already uploaded! It must be something wrong in

the code!



______________________________________________

R-help@r-project.org<mailto:R-help@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.






        [[alternative HTML version deleted]]



______________________________________________
R-help@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.

Reply via email to