Julia 中多个图像或数组的最大值

Question

我想找到最多几个图像：将它们加载到一个数组中并在第一维上找到一个最大值。

Python 代码示例：

import cv2
import sys
import numpy as np

imgs_paths = sys.argv[1:]
imgs = list(map(cv2.imread, imgs_paths))
imgs_arr = np.array(imgs, dtype=np.float32)
imgs_max = np.max(imgs_arr, 0)

我所做的如下：

using Colors, Images

function im_to_array(im)
    img_array = permutedims(channelview(im), (2,3,1)) 
    img_array = Float32.(img_array)
    return img_array
end


imgs = map(Images.load, imgs_paths)
imgs_arr = map(im_to_array, imgs)
a = imgs_arr
b = reshape(cat(a..., dims=1), tuple(length(a), size(a[1])...))
imgs_max = maximum(b, dims=1)

但它很丑。

我找到了更简单的方法来获得最大值（代码如下），但它的性能很糟糕。可能不是我所期待的。

function im_to_array(im)
    img_array = permutedims(channelview(im), (2,3,1)) 
    img_array = Float32.(img_array)
    return img_array
end

imgs = map(Images.load, imgs_paths)
imgs_arr = map(im_to_array, imgs)
imgs_max = max.(imgs_arr...)

第一种方法在 120 个全高清图像上的运行时间在我的笔记本电脑上约为 5 秒。而且我无法确定第二种方法的运行时间，因为我等待了大约 30 分钟并且它没有停止。我在 Julia 1.4.1 上测试它

有没有更好的方法来查找最多多个图像？

UPD：这是我想要的简单案例：

a = [zeros(Int8, 8, 8, 3), zeros(Int8, 8, 8, 3), zeros(Int8, 8, 8, 3)] # 3 black images with shape 8x8
max.(a) #doesn't work
max.(a...) #works with this simple input but when I test it on 120 FHD images it's extremely slow 

UPD2：我在少量图像上测试了这两种方法。

function max1(imgs_arr)
    a = imgs_arr
    b = reshape(cat(a..., dims=1), tuple(length(a), size(a[1])...))
    imgs_max = maximum(b, dims=1)
    return imgs_max
end

function max2(imgs_arr)
    return max.(imgs_arr...)
end

imgs_arr = my_imgs_arrays[1:5]

@time max1(imgs_arr)
@time max2(imgs_arr)

  0.247060 seconds (5.29 k allocations: 142.657 MiB)
  0.154158 seconds (44.85 k allocations: 26.388 MiB)

imgs_arr = my_imgs_arrays[1:15]

@time max1(imgs_arr)

@time max2(imgs_arr)

  0.600093 seconds (72.38 k allocations: 382.923 MiB)
  0.769446 seconds (1.24 M allocations: 71.374 MiB)

imgs_arr = my_imgs_arrays[1:25]

@time max1(imgs_arr)

@time max2(imgs_arr)

  1.057548 seconds (23.08 k allocations: 618.309 MiB)
  5.270050 seconds (151.52 M allocations: 2.329 GiB, 4.77% gc time)

所以，我使用的图片越多，效果就越慢。

Answer 1

您似乎希望对多个图像进行成对的最大缩减。首先，这是一个生成随机“图像”的函数：

rand_images(k, dims...) = [rand(UInt8, dims...) for _ = 1:k]

我将生成三个随机 10x12 图像的向量：

julia> images = rand_images(3, 10, 12)
3-element Array{Array{UInt8,2},1}:
 [0x51 0xdc … 0xf7 0x1e; 0xe1 0x10 … 0xd8 0x98; … ; 0x54 0x45 … 0x7a 0xaf; 0x7b 0xfc … 0x0a 0x81]
 [0xc8 0xa5 … 0xa8 0x81; 0x92 0x89 … 0x9f 0xbe; … ; 0x6a 0x03 … 0xb1 0xfd; 0x34 0xa9 … 0xa3 0x50]
 [0x26 0x9b … 0x2a 0x7c; 0x5c 0x7d … 0x8d 0x2b; … ; 0x32 0x1b … 0x57 0xdf; 0x96 0xa1 … 0x2a 0xc9]

一种直接的方法是进行成对的最大缩减：

julia> using BenchmarkTools

julia> @btime reduce(images) do a, b
           max.(a, b)
       end
  400.485 ns (2 allocations: 416 bytes)
10×12 Array{UInt8,2}:
 0xc8  0xdc  0x82  0xa7  0xa6  0xce  0xcd  0xb2  0x6e  0xba  0xf7  0x81
 0xe1  0x89  0x9f  0xeb  0x89  0xdf  0xd2  0xd2  0xab  0xea  0xd8  0xbe
 0xeb  0xdd  0x9e  0xe2  0xf5  0x4b  0xd2  0xe8  0xe4  0xf8  0xb9  0xf8
 0x63  0xa3  0xd7  0xea  0xf0  0x93  0xed  0xf7  0xfb  0xfb  0x9f  0xbb
 0xf2  0x51  0xf0  0xd4  0xfc  0xcf  0xf4  0xdd  0xeb  0xc3  0xe9  0xf9
 0xf8  0x72  0xfa  0x92  0x72  0xaa  0xa2  0xed  0xa1  0xdf  0xf1  0xd0
 0xef  0xe6  0x64  0xb3  0xd0  0x6a  0xce  0x9e  0x96  0xba  0xed  0xf9
 0xdb  0xc5  0x52  0xb3  0xf7  0xd1  0xdd  0xba  0xac  0xbc  0xd3  0xa1
 0x6a  0x45  0x88  0xda  0xf5  0xc6  0xcf  0x64  0xbc  0xf9  0xb1  0xfd
 0x96  0xfc  0xb1  0xc0  0xc4  0xcf  0x89  0xb4  0xe8  0xad  0xa3  0xc9

这相当快：400ns。我会在尺寸与你正在做的事情相当的图像上计时，但你没有提到我可以看到的图像尺寸（代码不依赖于数据，所以图像中的数据应该无关紧要）。

减少计算一个最大切片，一次减少一个图像，这可能不是最快的方法。似乎在所有图像中一次计算每个最大“像素”可能会更快，这有点复杂但也可以完成：

function max_images(images::Vector{<:Array})
    M = copy(images[1])
    for i = 1:length(M)
        for j = 2:length(images)
            M[i] = max(M[i], images[j][i])
        end
    end
    return M
end

这可行，但它需要 421 纳秒，这比数组缩减版本慢！哎呀。原因之一是无法保证图像的大小都相同，因此在每个图像的内部循环索引中都有边界检查。我们可以通过在@inbounds M[i] = max(M[i], images[j][i]) 上添加入站注释来跳过它，风险自负。这将时间降低到 282 ns。通过告诉编译器它可以安全地重新排序两个循环以通过在每个 for 循环上放置 @simd 宏来利用指令级并行性，可以获得更快的速度。这将时间缩短到 240 ns。最终版本的代码是：

function max_images(images::Vector{<:Array})
    M = copy(images[1])
    @simd for i = 1:length(M)
        @simd for j = 2:length(images)
            @inbounds M[i] = max(M[i], images[j][i])
        end
    end
    return M
end

Answer 2

我已经习惯了在 numpy 中避免循环，以至于我忘记了我可以使用它们。

我可以使用简单的循环（下面的代码中的max3），如@stefankarpinski 回答：

function max1(imgs_arr)
    a = imgs_arr
    b = reshape(cat(a..., dims=1), tuple(length(a), size(a[1])...))
    imgs_max = maximum(b, dims=1)
    return imgs_max
end

function max2(imgs_arr)
    return max.(imgs_arr...)
end

function max3(images::Vector{<:Array})
    M = copy(images[1])
    @simd for j = 2:length(images)
        M = max.(M, images[j])
    end
    return M
end

这是最快的方法：

# typeof(my_imgs) is Array{Array{Float32,3},1}
# size(my_imgs[1]) is (1080, 1920, 3)

imgs_arr = my_imgs[1:20] 
@time max1(imgs_arr)
@time max2(imgs_arr)
@time max3(imgs_arr)

  0.656771 seconds (5.62 k allocations: 498.630 MiB)
  3.237826 seconds (118.20 M allocations: 1.784 GiB, 5.24% gc time)
  0.137279 seconds (40 allocations: 474.611 MiB)

但是max 是简单的函数，我的问题仍然是其他函数，例如median。