GeoArrays produces positively y spaced rasters by default

[kongdd]

I am writing tiff using the following scripts:

begin
    ga = GeoArray(rand(360, 150))
    bbox!(ga, (min_x=-180.0, min_y=-60.0, max_x=180.0, max_y=90.0))
    epsg!(ga, 4326)  # in WGS84
    GeoArrays.write!("test.tif", ga)    
end

When reading in R, the crsTransform is error. The latitude range is changed to [-59, 91].

library(raster)
b <- brick("test.tif")

Warning message:
In .rasterFromGDAL(x, band = band, objecttype, ...) :
  data seems flipped. Consider using: flip(x, direction='y')

b
class      : RasterBrick 
dimensions : 150, 360, 54000, 1  (nrow, ncol, ncell, nlayers)
resolution : 1, 1  (x, y)
extent     : -180, 180, -59, 91  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : N:/Research/GEE_repos/GEE-latest/test.tif 
names      : layer 

[kongdd]

This problem can be solve by fixing bbox_to_affine. https://github.com/geo-julia/GeoArrays.jl/commit/2cc1c0912dd4d285bea3b7bf2135160a376b8a7a

[evetion]

Thanks for reporting! I agree that the data is flipped in the GDAL ecosystem sense, as it expects a negative y spacing.

But it shouldn't be wrong, especially not in terms of things like latitude:

❯ gdalinfo test.tif
Driver: GTiff/GeoTIFF
Files: test.tif
Size is 360, 150
Coordinate System is:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
Origin = (-180.000000000000000,-60.000000000000000)
Pixel Size = (1.000000000000000,1.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=DEFLATE
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (-180.0000000, -60.0000000) (180d 0' 0.00"W, 60d 0' 0.00"S)
Lower Left  (-180.0000000,  90.0000000) (180d 0' 0.00"W, 90d 0' 0.00"N)
Upper Right ( 180.0000000, -60.0000000) (180d 0' 0.00"E, 60d 0' 0.00"S)
Lower Right ( 180.0000000,  90.0000000) (180d 0' 0.00"E, 90d 0' 0.00"N)
Center      (   0.0000000,  15.0000000) (  0d 0' 0.01"E, 15d 0' 0.00"N)
Band 1 Block=256x256 Type=Float64, ColorInterp=Gray

It boils down to how you would interpret an Array in Julia. Is it positively spaced? For your fix to work, the data needs to be flipped upside down too! Note there's a GeoArrays.flipud! utility exactly for this. I see a few possibilities:

• We give a warning only when positive y spacing is encountered, hinting to use flipud!.
• Always use flipud!, to arrive at the negative y spacing.

@visr Any thoughts on this?

[visr]

GDAL seems to support both positive and negative y spacing, so I think the current behaviour is fine? Although negative y spacing seems to be the norm. Are there many applications that only work with negative cell spacing? Does QGIS show it well?

[evetion]

Showing is often no problem, but some GDAL command line tools refuse to work with them.