Rasters.jl icon indicating copy to clipboard operation
Rasters.jl copied to clipboard

`cellsize` should error informatively if there is no CRS

Open rafaqz opened this issue 1 year ago • 2 comments

This isn't so nice:

julia> lon = X(-179.875:1:179.875);

julia> lat = Y(-89.875:1:89.875);

julia> ras = Raster(rand(lon, lat));

julia> Rasters.cellsize(ras)
ERROR: MethodError: Cannot `convert` an object of type Nothing to an object of type CoordSys

Closest candidates are:
  convert(::Type{T1}, ::T2) where {T1<:GeoFormat, T2<:T1}
   @ GeoFormatTypes ~/.julia/packages/GeoFormatTypes/UFNTK/src/GeoFormatTypes.jl:93
  convert(::Type{T}, ::T) where T
   @ Base Base.jl:84
  convert(::Type{T}, ::AbstractString) where T<:GeoFormat
   @ GeoFormatTypes ~/.julia/packages/GeoFormatTypes/UFNTK/src/GeoFormatTypes.jl:154
  ...

Stacktrace:
 [1] cellsize(dims::Tuple{X{DimensionalData.Dimensions.Lookups.Sampled{…}}, Y{DimensionalData.Dimensions.Lookups.Sampled{…}}})
   @ RastersArchGDALExt ~/.julia/dev/Rasters/ext/RastersArchGDALExt/cellsize.jl:54
 [2] cellsize(x::Raster{Float64, 2, Tuple{X{…}, Y{…}}, Tuple{}, Matrix{Float64}, DimensionalData.NoName, DimensionalData.Dimensions.Lookups.NoMetadata, Nothing})
   @ RastersArchGDALExt ~/.julia/dev/Rasters/ext/RastersArchGDALExt/cellsize.jl:83
 [3] top-level scope
   @ REPL[70]:1
Some type information was truncated. Use `show(err)` to see complete types.

@tiemvanderdeure

rafaqz avatar May 08 '24 12:05 rafaqz

Surely there are a bunch of other functions that should behave similarly? They should probably all throw the same error? Or in some cases warn and assume WGS84?

For example this weird resample just gives a warning

lon = X(-179.875:1:179.875);
lat = Y(-89.875:1:89.875);
ras = Raster(rand(lon, lat));
resample(ras; crs = EPSG(4326))
┌ Warning: You have set a crs to resample to, but the object does not have crs so GDAL will assume it is already in the target crs. Use `newraster = setcrs(raster, somecrs)` to fix this.
└ @ RastersArchGDALExt C:\Users\tsh371\.julia\packages\Rasters\GI5GA\ext\RastersArchGDALExt\resample.jl:60

But reproject(ras; crs = EPSG(4326)) does not even throw a warning and just returns ras.

tiemvanderdeure avatar May 08 '24 13:05 tiemvanderdeure

Surely there are more! This is just one that came up for someone. Something consistent would be good.

rafaqz avatar May 08 '24 15:05 rafaqz