rrasternetcdfmap-projectionsreprojection-error

Reprojection oceancolor data in raster?


regarding my previous question Read/Open Oceancolor data from MODIS with nc format. I have tried to read the ncdf oceancolor data download here using ncdf4 package with this code

library(ncdf4)
nc <- nc_open('A2014325053500.L2_LAC_SST.nc')

# Get data for each variabel name as a large matrix
sst <- ncvar_get(nc, varid = 'geophysical_data/sst')
lon <- ncvar_get(nc, varid = 'navigation_data/longitude')
lat <- ncvar_get(nc, varid = 'navigation_data/latitude')

than this data, I convert to raster. But I failed to reproject the raster figure 1. The image should be like figure 2

library(sp)
library(raster)

r <- raster(sst,
        xmn=min(lat), xmx=max(lat), 
        ymn=min(lon), ymx=max(lon),
        crs=CRS('+proj=longlat +ellps=WGS84 +datum=WGS84'))

proj='+proj=longlat +ellps=WGS84 +datum=WGS84'
r2 <- projectRaster(r, crs=proj, method = 'ngb')
plot(r2)

Based on jbaums first suggestion the result is not correctly reprojected figure 3

I just realize that there is something wrong with the result when I increase the resolution. In the western and eastern part become white (no data). I used this code below:

r <- rasterize(df, raster(extent(df), res=0.0103011), 'sst', fun=mean)

Figure 1 figure 1

Figure 2 figure 2

Figure 3 Figure 3

Figure 4 Figure 4


Solution

  • In their current projection, the coordinates don't fall on a regular grid.

    One way around this is to represent the coordinates as a SpatialPoints object, and then, if you need a raster, rasterize it.

    library(ncdf4)
    library(raster)
    
    # Read in the data
    nc <- nc_open('~/../Downloads/A2014325053500.L2_LAC_SST.nc')
    sst <- ncvar_get(nc, varid = 'geophysical_data/sst')
    lon <- ncvar_get(nc, varid = 'navigation_data/longitude')
    lat <- ncvar_get(nc, varid = 'navigation_data/latitude')
    nc_close(nc)
    
    # Create a SpatialPointsDataFrame
    p <- data.frame(lon=c(lon), lat=c(lat), sst=c(sst))
    coordinates(p) <- ~lon+lat
    proj4string(p) <- '+init=epsg:4326'
    
    # Rasterize with appropriate resolution and aggregation function
    r <- rasterize(p, raster(extent(p) * 1.04, res=0.05), 'sst', fun=mean)
    # for high resolution look at ?gdalUtils::gdal_rasterize for efficiency
    

    Here's a plot:

    library(rasterVis)
    levelplot(r, at=seq(-75, 75, length=100), margin=FALSE, par.settings=BuRdTheme,
                  colorkey=list(height=0.6), main='MODIS Sea Surface Temperature')
    

    enter image description here