Plotting netCDF4 data with the right projection

Question

I guess this will sound simple to you but I just cannot get the right plot within R. I have a netcdf4 file and would like to visualize it within R just like I would do in QGIS and ArcGIS.

sample data

I tried to use insights from similar existing posts but did not succeed. The image get squeezed as shown below:

library(raster)
mycrs <- CRS("+proj=lcc +lat_1=35 +lat_2=51 +lat_0=39 +lon_0=14 +k=0.684241 +units=m +datum=WGS84 +no_defs")
pr <- raster(file.nc4, varname="precipitationCal")
projection(pr) <- mycrs
pr

pr

class       : RasterLayer 
dimensions  : 3600, 1800, 6480000  (nrow, ncol, ncell)
resolution  : 0.1, 0.1  (x, y)
extent      : -90, 90, -180, 180  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=lcc +lat_1=35 +lat_2=51 +lat_0=39 +lon_0=14 +k=0.684241 +units=m +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0 
data source : C:\Users\3B-HHR.MS.MRG.3IMERG.20150405-S170000-E172959.1020.V03D.HDF5.dap.nc4 
names       : precipitationCal 
zvar        : precipitationCal 

plot(pr)

NB: my region of interest is lat=c(41,84) and lon=c(-142,-53)

file.nc4 attributes:

Grid.GridHeader: BinMethod=ARITHMETIC_MEAN;
Registration=CENTER;
LatitudeResolution=0.1;
LongitudeResolution=0.1;
NorthBoundingCoordinate=90;
SouthBoundingCoordinate=-90;
EastBoundingCoordinate=180;
WestBoundingCoordinate=-180;
Origin=SOUTHWEST;

EDIT: Using @mdsummer's suggestion below:

r <- '3B-HHR.MS.MRG.3IMERG.20150405-S170000-E172959.1020.V03D.HDF5.dap.nc4'
pr <- raster(r, varname="precipitationCal")
library(maptools);
data(wrld_simpl)
projection(wrld_simpl)
proj4string(pr)=proj4string(wrld_simpl)
plot(t(pr))
plot(wrld_simpl, add = TRUE)

I got this image but it is still not on the right latlon projection?.

Answer 1

The problem is you are assigning metadata to the object, and probably even overwriting the existing metadata. This happens here:

projection(pr) <- mycrs

You should first run this to see what raster thinks it is already:

projection(pr)

or just print the object out to get a fully summary:

pr

Going out on a limb, I think you should warp this raster to get to your target:

library(raster)
mycrs <- CRS("+proj=lcc +lat_1=35 +lat_2=51 +lat_0=39 +lon_0=14 +k=0.684241 +units=m +datum=WGS84 +no_defs")
pr <- raster(file.nc4, varname="precipitationCal")

ppr <- projectRaster(pr, crs = mycrs)  ## create a new object for good measure

BUT, warping is a destructive no-holds-barred re-modelling process, and so you may need to limit the extent of the target or crop the input or many other possible options. It also looks like your data is transpose y/x and so you probably need to fix that first - make sure you can get a sensible map, try this:

 pr <- raster(file.nc4, varname="precipitationCal")
 plot(t(pr))
 library(maptools);data(wrld_simpl)
 plot(wrld_simpl, add = TRUE)

(Untested since I can't access dropbox files right now).

If raster didn't already know what projection your data is in ("+proj=longlat +ellps=WGS84 etc. etc") then the projectRaster step would fail, and you would need to assign the right longlat PROJ.4 string before reprojecting.