DimensionMismatch for a Raster

Specific Domains Geo
1

A Raster behaves just like an Array in suitable contexts, but I found the following difference in broadcast behavior. Is it 1) an intended behavior, 2) an unavoidable inconsistency, or 3) a bug?

Basically,

avr = sum(v,dims=3)/size(v,3)
dev = v .- avr

works when v is a regular 3D array but it doesn’t when v is a Raster. With a Raster, you need avr[:,:,1]. Below is a self-contained example.

I’m using Julia 1.8.0 on macOS 12.5.1 .

using Downloads
using Rasters

const url = "https://www.unidata.ucar.edu/software/netcdf/examples/sresa1b_ncar_ccsm3-example.nc"
infile = Downloads.download(url, "tmp.nc")
#infile = "tmp.nc"

ua = Raster(infile; name = "ua")

v = ua[:,:,:,1]
avr = sum(v,dims=3)/size(v,3)
#dev = v .- avr # => size mismatch!
dev = v .- avr[:,:,1] # => Fine.

v = Array{Union{Missing,Float64}}(undef, size(ua[:,:,:,1]))

v .= ua[:,:,:,1]
avr = sum(v,dims=3)/size(v,3)
dev = v .- avr
2

This is a bug in DimensionalData.jl. But a very subtle one.

For broadcast we need to compare that the dimensions are the same (DD promises not to broadcast different dimensions together or the returned dimensions are ambiguous), but we need to relax the comparison a little.

Basically we will need to allow mismatches where the length of one dimension is 1, and make sure to keep the longer dimension.

1 Like