Calculate the average of 2 'image's in R when the resolution is not the same

Question

I have an image in R that is formed by X and Y plus a Z value which is in a matrix format. I use image to give me the output of that. I'd like to calculate the average value of the 2 images so I can join them together. There might be a way to use spatial packages to do that, but can't figure out a way to do this. The major problem is that the 2 Z matrices are not having the same number of columns, making the comparison between then difficult.

par(mfrow=c(1,1))
x1 = seq(1,10, by =1)
y1 = seq(1,10, by =1)
z1 = outer(x1,y1)

x2 = seq(7,12, by =.01)
y2 = seq(5,12, by =.01)
z2 = outer(x2,y2, FUN = "*")

image(x1,y1,z1, xlim = range(c(x1,x2)), ylim = range(c(y1,y2)), asp = 1)
image(x2,y2,z2, add =TRUE)

get.x.1 = x1 > min(x2)
get.y.1 = y1 > min(y2)
get.x.2 = x2 < max(x1)
get.y.2 = y2 < max(y1)

segments(min(x2),min(y2),max(x1),min(y2), lwd = 2)
segments(max(x1),min(y2),max(x1),max(y1), lwd = 2)
segments(max(x1),max(y1),min(x2),max(y1), lwd = 2)
segments(min(x2),max(y1),min(x2),min(y2), lwd = 2)

I'd like to average the 2 images below where they overlap (in the black rectangle that was added using the segments)

1. first without a buffer and
2. with a buffer (so that the buffer takes into account more surrounding values as a way to smooth the averaging process).

The end result should be a combined Z matrix with the average part but also the parts of the image that is not overlapping.

Answer 1

It appears that magick is not agnostic on the importance of 'dim'[s], albeit more compact:

library(magick)
Linking to ImageMagick 7.1.0.46
Enabled features: cairo, fontconfig, freetype, fftw, heic, lcms, pango, raw, rsvg, webp, x11
Disabled features: ghostscript
Using 4 threads
# again mangling your z1 z2
z1 <- z1/100
z1_rast <- as.raster(z1) # here grDevices::as.raster
z1_mgk <- image_read(z1_rast)
z2_1k <- z2/1000
z2_1k <- z2_1k + 0.5
z2_1k_rast <- as.raster(z2_1k)
z2_mgk <- image_read(z2_1k_rast)
happy_z12_2 <- image_composite(z2_mgk, z1_mgk, operator = 'blend', offset='+500+x150', compose_args = '50')
plot(happy_z12_2)
# a discerning eye can just make out the 10x10 'atop'and trust to 50%

Someone more experienced with magick will likely present a better approach. I was happy to get ImageMagic-7 built.

Answer 2

Perhaps a terra approach with glaring deviations from your presented data that might not serve your workflow, and in this case will result in an approximation of 'buffered'. In either case below, the difference in dim must be addressed:

library(terra)

x1a = seq(1,10, by =.01) # first glaring deviation
y1a = seq(1,10, by =.01) # and second
z1a = outer(x1a,y1a)

x1b = seq(1,10,by = 1)
y1b = seq(1,10,by = 1)
z1b = outer(x1b,y1b)

x2 = seq(7,12, by =.01)
y2 = seq(5,12, by =.01)
z2 = outer(x2,y2, FUN = "*")

z1a_apply <- apply(z1a, 2, FUN = 'rev') # get value gradients reversed
z2_apply <- apply(z2, 2, FUN = 'rev') # get value gradients reversed

z1a_rast <- rast(z1a_apply)
z2_rast <- rast(z2_apply)
# these leave origins at (0,0) which is not the case

ext(z2_rast) <- c(701, 1402, 501, 1002)
z1a_z2_crop_ext <- ext(crop(z1a_rast, z2_rast))
z1a_crop <- crop(z1a_rast, z1a_z2_crop_ext)
z2_crop <- crop(z2_rast, z1a_z2_crop_ext)
z1a_z2_mean <- app(c(z1a_crop, z2_crop), mean)
z_sprc <- sprc(z1a_z2_mean, z1a_rast, z2_rast)
z_merge <- merge(z_sprc)
plot(z_merge)

And as will be seen, my fontconfig is broken. So, a partial approach so far.

# picking up with z1b
x1b = seq(1,10,by = 1)
y1b = seq(1,10,by = 1)
z1b = outer(x1b,y1b)
z1b_rast <- rast(z1b)
dim(z1b_rast) <- c(1000, 1000)
# z1b_rast has lost all values
values(z1b_rast) <- outer(c(1:1000),c(1:1000))
# but here the gradient is wrong with high values lower right
z1b_flip <- flip(z1b_rast)
# a picture of a cat might not survive this treatment
# extent, resolution, and origin also have to be adjusted
ext(z1b_flip) <- c(0,1000,0,1000)
res(z1b_flip) <- c(1,1)
origin(z1b_flip) <- c(0.5, 0.5)# should have been done on z2

The subtleties in the 'polygon' approach via line seqments to be addressed next. To better approximate your plot above, points are inset from your segments above by -50.

library(sf)
# using your segments
pts1 <- matrix(c(min(x2),min(y2),max(x1a),min(y2)), nrow = 2, byrow = TRUE)*100
pts2 <- matrix(c(max(x1a),min(y2),max(x1a),max(y1a)), nrow = 2, byrow = TRUE)*100
pts3 <- matrix(c(max(x1a),max(y1a),min(x2),max(y1a)),nrow =2, byrow = TRUE)*100
pts4 <- matrix(c(min(x2),max(y1a),min(x2),min(y2)),nrow = 2, byrow=TRUE)*100
# one point from each line, all inset by 50
pts1_2 <- pts1[2, ] -50
pts2_2 <- pts2[2, ] -50
pts3_2 <- pts3[2, ] -50
pts4_2 <- pts4[2, ] -50

We later find this isn't actually useful compensation as the resulting poly falls outsize our z2_rast, and to avoid general madness it is better to have a poly that is within both rast(s) to pull our mean.

And so adjust out points...

pts1_2[2] <- pts1_2[2] + 75
pts4_2 <- pts4_2 + 75
pts3_2[1] <- pts3_2[1] + 75
poly4 <- st_cast(c(st_point(pts1_2), st_point(pts2_2), st_point(pts3_2), st_point(pts4_2)), 'POLYGON')
#make poly4 a SpatVector
poly4_vect <- vect(poly4)
z1b_poly4_crop <- crop(z1b_flip, poly4_vect)
z2_poly4_crop <- crop(z2_rast, poly4_vect)
#check for damages
all.equal(dim(z1b_poly4_crop), dim(z2_poly4_crop))
[1] TRUE
# Oh, hurray!, But
z1b_z2_poly4_mean <- app(c(z1b_poly4_crop, z2_poly4_crop), mean)
Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 'app': [rast] extents do not match
ext(z1b_poly4_crop)
SpatExtent : 725.5, 950.5, 525.5, 950.5 (xmin, xmax, ymin, ymax)

ext(z2_poly4_crop)
SpatExtent : 725, 950, 525, 950 (xmin, xmax, ymin, ymax)
# so adjust origin on z2_poly4_crop, or better, on z2 before...
origin(z2_rast)
[1] 0 0
origin(z2_rast) <- c(0.5, 0.5)
z2_poly4_crop <- crop(z2_rast, poly4_vect)

ext(z2_poly4_crop)
SpatExtent : 725.5, 950.5, 525.5, 950.5 (xmin, xmax, ymin, ymax)
# and now can pull mean
z1b_z2_poly4_mean <- app(c(z1b_poly4_crop, z2_poly4_crop), mean)
###
mean_z1b_flip_z2 <- sprc(z1b_z2_poly4_mean, z1b_flip, z2_rast)
mean_z1b_flip_z2_sprc <- sprc(z1b_z2_poly4_mean, z1b_flip, z2_rast)
mean_flip_z2_merge <- merge(mean_z1b_flip_z2_sprc)
# or mosaic - mean goes last and plot

Lots of good fun and rabbit holes to stub one's toes on. There are likely much more compact approaches that others might offer. I imagine that much of this could also be approached through magick.