converting raster to netcdf .nc with keeping CRS information

classic Classic list List threaded Threaded
5 messages Options
Reply | Threaded
Open this post in threaded view
|

converting raster to netcdf .nc with keeping CRS information

Ahmed El-Gabbas

Hello all,


I am having a problem while converting a raster object as NetCDF (.nc) file, with keeping the CRS information in the output file.

Here is a reproducible code:

require(raster)
require(ncdf4)
CurrTemp <- tempfile()
download.file(url = "https://seaice.uni-bremen.de/data/amsre/asi_daygrid_swath/s6250/2003/feb/Antarctic/asi-s6250-20030214-v5.hdf", destfile = CurrTemp, mode = "wb", quiet = T)
r <- raster(CurrTemp)
r <- flip(r,2)
extent(r) <- c(-3950000, 3950000, -3950000, 4350000)
crs(r) <- "+proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs"
r
# class       : RasterLayer 
# dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
# resolution  : 6250, 6250  (x, y)
# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs 
# data source : in memory
# names       : layer 
# values      : 0, 100  (min, max)

So far the raster object reads well, with CRS information included.

However, when I try to save it as .nc file, R prints "coord. ref. : NA", and the produced file does not contain the CRS information.

writeRaster(r, filename = "O:/Ahmed/r001.nc", varname="IceConc",
            overwrite=TRUE, format="CDF",
            xname="Longitude", yname="Latitude")
# class       : RasterLayer 
# dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
# resolution  : 6250, 6250  (x, y)
# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
# coord. ref. : NA 
# data source : O:/Ahmed/r001.nc 
# names       : IceConc 
# zvar        : IceConc

raster("O:/Ahmed/r001.nc")
# class       : RasterLayer 
# dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
# resolution  : 6250, 6250  (x, y)
# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
# coord. ref. : NA 
# data source : O:/Ahmed/r001.nc 
# names       : IceConc 
# zvar        : IceConc

Any solution?


N.B. I also sent the same question at stackoverflow, without a solution so far.


Thanks


~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Dr. Ahmed El-Gabbas,
Ocean Acoustics Lab,
Alfred-Wegener-Institut
Helmholtz-Zentrum für Polar und Meeresforschung
image.png

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting raster to netcdf .nc with keeping CRS information

Ben Tupper
Hi,

Have you tried using the 'prj' argument to writeRaster?  I don't know that it will be the solution, but according to ...

https://www.rdocumentation.org/packages/raster/versions/2.8-19/topics/writeRaster <https://www.rdocumentation.org/packages/raster/versions/2.8-19/topics/writeRaster>

It is documented among the '...' arguments.  Setting prj = TRUE will cause the crs to be written to the file.

Cheers,
Ben

> On Feb 12, 2019, at 8:34 AM, Ahmed El-Gabbas <[hidden email]> wrote:
>
> Hello all,
>
> I am having a problem while converting a raster object as NetCDF (.nc) file, with keeping the CRS information in the output file.
> Here is a reproducible code:
> require(raster)
> require(ncdf4)
> CurrTemp <- tempfile()
> download.file(url = "https://seaice.uni-bremen.de/data/amsre/asi_daygrid_swath/s6250/2003/feb/Antarctic/asi-s6250-20030214-v5.hdf <https://seaice.uni-bremen.de/data/amsre/asi_daygrid_swath/s6250/2003/feb/Antarctic/asi-s6250-20030214-v5.hdf>", destfile = CurrTemp, mode = "wb", quiet = T)
> r <- raster(CurrTemp)
> r <- flip(r,2)
> extent(r) <- c(-3950000, 3950000, -3950000, 4350000)
> crs(r) <- "+proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs"
> r
> # class       : RasterLayer
> # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
> # resolution  : 6250, 6250  (x, y)
> # extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
> # coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs
> # data source : in memory
> # names       : layer
> # values      : 0, 100  (min, max)
> So far the raster object reads well, with CRS information included.
> However, when I try to save it as .nc file, R prints "coord. ref. : NA", and the produced file does not contain the CRS information.
> writeRaster(r, filename = "O:/Ahmed/r001.nc <http://r001.nc/>", varname="IceConc",
>             overwrite=TRUE, format="CDF",
>             xname="Longitude", yname="Latitude")
> # class       : RasterLayer
> # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
> # resolution  : 6250, 6250  (x, y)
> # extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
> # coord. ref. : NA
> # data source : O:/Ahmed/r001.nc <http://r001.nc/>
> # names       : IceConc
> # zvar        : IceConc
>
> raster("O:/Ahmed/r001.nc <http://r001.nc/>")
> # class       : RasterLayer
> # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
> # resolution  : 6250, 6250  (x, y)
> # extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
> # coord. ref. : NA
> # data source : O:/Ahmed/r001.nc <http://r001.nc/>
> # names       : IceConc
> # zvar        : IceConc
> Any solution?
>
> N.B. I also sent the same question at stackoverflow <https://stackoverflow.com/questions/54593552/saving-r-raster-to-netcdf-nc-file-with-keeping-crs-information>, without a solution so far.
>
> Thanks
>
> ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> Dr. Ahmed El-Gabbas,
> Ocean Acoustics Lab,
> Alfred-Wegener-Institut
> Helmholtz-Zentrum für Polar und Meeresforschung
> <image.png>
> My Website:  https://elgabbas.github.io <https://elgabbas.github.io/>_______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Ben Tupper
Bigelow Laboratory for Ocean Sciences
60 Bigelow Drive, P.O. Box 380
East Boothbay, Maine 04544
http://www.bigelow.org

Ecological Forecasting: https://eco.bigelow.org/






        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting raster to netcdf .nc with keeping CRS information

Ahmed El-Gabbas
Thanks Ben for your suggestion.

Adding *prj = TRUE *argument caused this error:
           Error in ncdf4::ncvar_def(varname, varunit, list(xdim, ydim),
NAvalue(x),  :
           unused argument (prj = T)

Ahmed

On Tue, Feb 12, 2019 at 3:13 PM Ben Tupper <[hidden email]> wrote:

> Hi,
>
> Have you tried using the 'prj' argument to writeRaster?  I don't know that
> it will be the solution, but according to ...
>
>
> https://www.rdocumentation.org/packages/raster/versions/2.8-19/topics/writeRaster
>
> It is documented among the '...' arguments.  Setting prj = TRUE will cause
> the crs to be written to the file.
>
> Cheers,
> Ben
>
> On Feb 12, 2019, at 8:34 AM, Ahmed El-Gabbas <[hidden email]> wrote:
>
> Hello all,
>
> I am having a problem while converting a raster object as NetCDF (.nc)
> file, with keeping the CRS information in the output file.
> Here is a reproducible code:
>
> require(raster)
> require(ncdf4)
> CurrTemp <- tempfile()
> download.file(url = "https://seaice.uni-bremen.de/data/amsre/asi_daygrid_swath/s6250/2003/feb/Antarctic/asi-s6250-20030214-v5.hdf", destfile = CurrTemp, mode = "wb", quiet = T)
> r <- raster(CurrTemp)
> r <- flip(r,2)
> extent(r) <- c(-3950000, 3950000, -3950000, 4350000)
> crs(r) <- "+proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs"
> r# class       : RasterLayer # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs # data source : in memory# names       : layer # values      : 0, 100  (min, max)
>
> So far the raster object reads well, with CRS information included.
> However, when I try to save it as .nc file, R prints "coord. ref. : NA",
> and the produced file does not contain the CRS information.
>
> writeRaster(r, filename = "O:/Ahmed/r001.nc", varname="IceConc",
>             overwrite=TRUE, format="CDF",
>             xname="Longitude", yname="Latitude")# class       : RasterLayer # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref. : NA # data source : O:/Ahmed/r001.nc # names       : IceConc # zvar        : IceConc
>
> raster("O:/Ahmed/r001.nc")# class       : RasterLayer # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref. : NA # data source : O:/Ahmed/r001.nc # names       : IceConc # zvar        : IceConc
>
> Any solution?
>
> N.B. I also sent the same question at stackoverflow
> <https://stackoverflow.com/questions/54593552/saving-r-raster-to-netcdf-nc-file-with-keeping-crs-information>,
> without a solution so far.
>
> Thanks
>
> ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> *Dr. Ahmed El-Gabbas,*
> Ocean Acoustics Lab,
> Alfred-Wegener-Institut
> Helmholtz-Zentrum für Polar und Meeresforschung
> <image.png>
> My Website:  https://elgabbas.github.io
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>
> Ben Tupper
> Bigelow Laboratory for Ocean Sciences
> 60 Bigelow Drive, P.O. Box 380
> East Boothbay, Maine 04544
> http://www.bigelow.org
>
> Ecological Forecasting: https://eco.bigelow.org/
>
>
>
>
>
>

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting raster to netcdf .nc with keeping CRS information

Ben Tupper
Hi Ahmed,

I can confirm that on my system I get the same result as you.

Does your workflow confine you to using NetCDF?  I can successfully write/read to the native raster format.

r2 <- writeRaster(r, filename = "r001.grd", overwrite=TRUE)
r2
#class       : RasterLayer
#dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
#resolution  : 6250, 6250  (x, y)
#extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs
#data source : /home/btupper/r001.grd
#names       : layer
#values      : 0, 100  (min, max)


r3 <- raster("r001.grd")
r3
#class       : RasterLayer
#dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
#resolution  : 6250, 6250  (x, y)
#extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs
#data source : /home/btupper/r001.grd
#names       : layer
#values      : 0, 100  (min, max)

Ben

> sessionInfo()
R version 3.5.1 (2018-07-02)
Platform: x86_64-redhat-linux-gnu (64-bit)
Running under: CentOS Linux 7 (Core)

Matrix products: default
BLAS/LAPACK: /usr/lib64/R/lib/libRblas.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8  
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C      

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base    

other attached packages:
[1] raster_2.8-19 sp_1.3-1      ncdf4_1.16  

loaded via a namespace (and not attached):
[1] compiler_3.5.1   rgdal_1.3-6      Rcpp_1.0.0       codetools_0.2-15
[5] grid_3.5.1       lattice_0.20-35





> On Feb 12, 2019, at 9:19 AM, Ahmed El-Gabbas <[hidden email]> wrote:
>
> Thanks Ben for your suggestion.
>
> Adding *prj = TRUE *argument caused this error:
>           Error in ncdf4::ncvar_def(varname, varunit, list(xdim, ydim),
> NAvalue(x),  :
>           unused argument (prj = T)
>
> Ahmed
>
> On Tue, Feb 12, 2019 at 3:13 PM Ben Tupper <[hidden email]> wrote:
>
>> Hi,
>>
>> Have you tried using the 'prj' argument to writeRaster?  I don't know that
>> it will be the solution, but according to ...
>>
>>
>> https://www.rdocumentation.org/packages/raster/versions/2.8-19/topics/writeRaster
>>
>> It is documented among the '...' arguments.  Setting prj = TRUE will cause
>> the crs to be written to the file.
>>
>> Cheers,
>> Ben
>>
>> On Feb 12, 2019, at 8:34 AM, Ahmed El-Gabbas <[hidden email]> wrote:
>>
>> Hello all,
>>
>> I am having a problem while converting a raster object as NetCDF (.nc)
>> file, with keeping the CRS information in the output file.
>> Here is a reproducible code:
>>
>> require(raster)
>> require(ncdf4)
>> CurrTemp <- tempfile()
>> download.file(url = "https://seaice.uni-bremen.de/data/amsre/asi_daygrid_swath/s6250/2003/feb/Antarctic/asi-s6250-20030214-v5.hdf", destfile = CurrTemp, mode = "wb", quiet = T)
>> r <- raster(CurrTemp)
>> r <- flip(r,2)
>> extent(r) <- c(-3950000, 3950000, -3950000, 4350000)
>> crs(r) <- "+proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs"
>> r# class       : RasterLayer # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs # data source : in memory# names     : layer # values      : 0, 100  (min, max)
>>
>> So far the raster object reads well, with CRS information included.
>> However, when I try to save it as .nc file, R prints "coord. ref. : NA",
>> and the produced file does not contain the CRS information.
>>
>> writeRaster(r, filename = "O:/Ahmed/r001.nc", varname="IceConc",
>>            overwrite=TRUE, format="CDF",
>>            xname="Longitude", yname="Latitude")# class       : RasterLayer # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref. : NA # data source : O:/Ahmed/r001.nc # names       : IceConc # zvar        : IceConc
>>
>> raster("O:/Ahmed/r001.nc")# class       : RasterLayer # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)# extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref. : NA # data source : O:/Ahmed/r001.nc # names       : IceConc # zvar        : IceConc
>>
>> Any solution?
>>
>> N.B. I also sent the same question at stackoverflow
>> <https://stackoverflow.com/questions/54593552/saving-r-raster-to-netcdf-nc-file-with-keeping-crs-information>,
>> without a solution so far.
>>
>> Thanks
>>
>> ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
>> *Dr. Ahmed El-Gabbas,*
>> Ocean Acoustics Lab,
>> Alfred-Wegener-Institut
>> Helmholtz-Zentrum für Polar und Meeresforschung
>> <image.png>
>> My Website:  https://elgabbas.github.io
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>>
>> Ben Tupper
>> Bigelow Laboratory for Ocean Sciences
>> 60 Bigelow Drive, P.O. Box 380
>> East Boothbay, Maine 04544
>> http://www.bigelow.org
>>
>> Ecological Forecasting: https://eco.bigelow.org/
>>
>>
>>
>>
>>
>>
>
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Ben Tupper
Bigelow Laboratory for Ocean Sciences
60 Bigelow Drive, P.O. Box 380
East Boothbay, Maine 04544
http://www.bigelow.org

Ecological Forecasting: https://eco.bigelow.org/

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: converting raster to netcdf .nc with keeping CRS information

Michael Sumner-2
I agree, I wouldn't try to bend NetCDF to work this way via raster. The
convention/s require parameters of the projection to be specified via
attributes and you could do that via RNetCDF or ncdf4, but it's a lot of
painful work.

It's not easy to create NetCDF from R, it's a simple as that - and
typically those files don't play well with projections - it's changed
recently but you'll see a lot of patchy partial support across different
softwares.

You can use rgdal or call out out to GDAL - but then round-tripping gives
warnings, because again raster has to unpack all the components of metadata
that GDAL created, and none of this is very systematic and clearly is
broken atm.

system('gdal_translate r001.nc r001_aux.nc -of NetCDF -a_srs
"+init=epsg:3412"')
## you could do this with  stars too
rgdal::writeGDAL(as(r, "SpatialGridDataFrame"), "r001_rgdal.nc", drivername
= "NetCDF")


But, beware as in neither case is raster guaranteed to get the right
interpretation (THIS IS WRONG):

raster("r001_aux.nc")
raster("r001_rgdal.nc")
class       : RasterLayer
dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
resolution  : 6250, 6250  (x, y)
extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=stere +lon_0=0 +x_0=0 +y_0=0 +lat_0=-90 +pm=0
+a=6378273 +rf=298.279411123064
data source : r001_rgdal.nc
names       : GDAL.Band.Number.1
zvar        : Band1

BTW, you called your X and Y "longitude" and "latitude" but this grid isn't
defined by those (eastings and northings is not a better name either, but
that's another problem).

Do you have to create new files?  I find that constantly a fraught problem,
For usage it's way easier to wrap stuff up with R code and work around all
these problems.

Cheers, Mike.

On Wed, 13 Feb 2019 at 12:47 Ben Tupper <[hidden email]> wrote:

> Hi Ahmed,
>
> I can confirm that on my system I get the same result as you.
>
> Does your workflow confine you to using NetCDF?  I can successfully
> write/read to the native raster format.
>
> r2 <- writeRaster(r, filename = "r001.grd", overwrite=TRUE)
> r2
> #class       : RasterLayer
> #dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
> #resolution  : 6250, 6250  (x, y)
> #extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin,
> ymax)
> #coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0
> +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs
> #data source : /home/btupper/r001.grd
> #names       : layer
> #values      : 0, 100  (min, max)
>
>
> r3 <- raster("r001.grd")
> r3
> #class       : RasterLayer
> #dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)
> #resolution  : 6250, 6250  (x, y)
> #extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin,
> ymax)
> #coord. ref. : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0
> +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs
> #data source : /home/btupper/r001.grd
> #names       : layer
> #values      : 0, 100  (min, max)
>
> Ben
>
> > sessionInfo()
> R version 3.5.1 (2018-07-02)
> Platform: x86_64-redhat-linux-gnu (64-bit)
> Running under: CentOS Linux 7 (Core)
>
> Matrix products: default
> BLAS/LAPACK: /usr/lib64/R/lib/libRblas.so
>
> locale:
>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
>  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8
>  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
>  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
>  [9] LC_ADDRESS=C               LC_TELEPHONE=C
> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
>
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods   base
>
> other attached packages:
> [1] raster_2.8-19 sp_1.3-1      ncdf4_1.16
>
> loaded via a namespace (and not attached):
> [1] compiler_3.5.1   rgdal_1.3-6      Rcpp_1.0.0       codetools_0.2-15
> [5] grid_3.5.1       lattice_0.20-35
>
>
>
>
>
> > On Feb 12, 2019, at 9:19 AM, Ahmed El-Gabbas <[hidden email]> wrote:
> >
> > Thanks Ben for your suggestion.
> >
> > Adding *prj = TRUE *argument caused this error:
> >           Error in ncdf4::ncvar_def(varname, varunit, list(xdim, ydim),
> > NAvalue(x),  :
> >           unused argument (prj = T)
> >
> > Ahmed
> >
> > On Tue, Feb 12, 2019 at 3:13 PM Ben Tupper <[hidden email]> wrote:
> >
> >> Hi,
> >>
> >> Have you tried using the 'prj' argument to writeRaster?  I don't know
> that
> >> it will be the solution, but according to ...
> >>
> >>
> >>
> https://www.rdocumentation.org/packages/raster/versions/2.8-19/topics/writeRaster
> >>
> >> It is documented among the '...' arguments.  Setting prj = TRUE will
> cause
> >> the crs to be written to the file.
> >>
> >> Cheers,
> >> Ben
> >>
> >> On Feb 12, 2019, at 8:34 AM, Ahmed El-Gabbas <[hidden email]>
> wrote:
> >>
> >> Hello all,
> >>
> >> I am having a problem while converting a raster object as NetCDF (.nc)
> >> file, with keeping the CRS information in the output file.
> >> Here is a reproducible code:
> >>
> >> require(raster)
> >> require(ncdf4)
> >> CurrTemp <- tempfile()
> >> download.file(url = "
> https://seaice.uni-bremen.de/data/amsre/asi_daygrid_swath/s6250/2003/feb/Antarctic/asi-s6250-20030214-v5.hdf",
> destfile = CurrTemp, mode = "wb", quiet = T)
> >> r <- raster(CurrTemp)
> >> r <- flip(r,2)
> >> extent(r) <- c(-3950000, 3950000, -3950000, 4350000)
> >> crs(r) <- "+proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0
> +y_0=0 +a=6378273 +b=6356889.449 +units=m +no_defs"
> >> r# class       : RasterLayer # dimensions  : 1328, 1264, 1678592
> (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)# extent      :
> -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref.
> : +proj=stere +lat_0=-90 +lat_ts=-70 +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378273
> +b=6356889.449 +units=m +no_defs # data source : in memory# names     :
> layer # values      : 0, 100  (min, max)
> >>
> >> So far the raster object reads well, with CRS information included.
> >> However, when I try to save it as .nc file, R prints "coord. ref. : NA",
> >> and the produced file does not contain the CRS information.
> >>
> >> writeRaster(r, filename = "O:/Ahmed/r001.nc", varname="IceConc",
> >>            overwrite=TRUE, format="CDF",
> >>            xname="Longitude", yname="Latitude")# class       :
> RasterLayer # dimensions  : 1328, 1264, 1678592  (nrow, ncol, ncell)#
> resolution  : 6250, 6250  (x, y)# extent      : -3950000, 3950000,
> -3950000, 4350000  (xmin, xmax, ymin, ymax)# coord. ref. : NA # data source
> : O:/Ahmed/r001.nc # names       : IceConc # zvar        : IceConc
> >>
> >> raster("O:/Ahmed/r001.nc")# class       : RasterLayer # dimensions  :
> 1328, 1264, 1678592  (nrow, ncol, ncell)# resolution  : 6250, 6250  (x, y)#
> extent      : -3950000, 3950000, -3950000, 4350000  (xmin, xmax, ymin,
> ymax)# coord. ref. : NA # data source : O:/Ahmed/r001.nc # names       :
> IceConc # zvar        : IceConc
> >>
> >> Any solution?
> >>
> >> N.B. I also sent the same question at stackoverflow
> >> <
> https://stackoverflow.com/questions/54593552/saving-r-raster-to-netcdf-nc-file-with-keeping-crs-information
> >,
> >> without a solution so far.
> >>
> >> Thanks
> >>
> >> ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
> >> *Dr. Ahmed El-Gabbas,*
> >> Ocean Acoustics Lab,
> >> Alfred-Wegener-Institut
> >> Helmholtz-Zentrum für Polar und Meeresforschung
> >> <image.png>
> >> My Website:  https://elgabbas.github.io
> >> _______________________________________________
> >> R-sig-Geo mailing list
> >> [hidden email]
> >> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >>
> >>
> >> Ben Tupper
> >> Bigelow Laboratory for Ocean Sciences
> >> 60 Bigelow Drive, P.O. Box 380
> >> East Boothbay, Maine 04544
> >> http://www.bigelow.org
> >>
> >> Ecological Forecasting: https://eco.bigelow.org/
> >>
> >>
> >>
> >>
> >>
> >>
> >
> >       [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > [hidden email]
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
> Ben Tupper
> Bigelow Laboratory for Ocean Sciences
> 60 Bigelow Drive, P.O. Box 380
> East Boothbay, Maine 04544
> http://www.bigelow.org
>
> Ecological Forecasting: https://eco.bigelow.org/
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
--
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo