Skip to main content
4 of 4
edited tags; edited title
Spacedman
  • 66.5k
  • 6
  • 83
  • 121

Raster file statistics in R (stars and raster package) don't match gdalinfo

I'm trying to read a gridded binary file in R using read_stars. Comparing the stats from stars package to the stats from gdalinfo, I could see differences in min/max values. Below is the sample code.

library(stars)
library(sf)
library(tidyverse)
library(fs)
library(raster)

url <- 'ftp://ftp.wpc.ncep.noaa.gov/pwpf/conus_2.5km/2.5kmpwpf_72hr/2.5kmprcntil_72hsnow_90pt_latestf072.grb'
file <- 'snow.grb'
download.file(url, destfile='snow.grb')

snow.stars <- read_stars('snow.grb')
snow.stars %>% as_tibble() -> snow.stars.tibble
snow.stars.tibble %>% summary(mean=mean(snow.grb, na.rm=T), min=min(snow.grb, na.rm=T), max=max(snow.grb, na.rm=T), sd=sd(snow.grb, na.rm=T))    

# snow.raster = raster('snow.grb')
# cellStats(snow.raster, stat='min', na.rm=T)
# cellStats(snow.raster, stat='max', na.rm=T)

with gdalinfo the stats are (gdalinfo -mm -stats snow.grb) Minimum=-0.000, Maximum=0.910, Mean=0.043, StdDev=0.115

The min/max values returned by gdalinfo are matching with qgis as well.

with stars, I'm getting these. max: 0.917, min: -1.17e-6, sd: 0.106 and mean: 0.0384

What am I missing?

Surya
  • 123
  • 4