对齐数组值

Question

假设我有两个数组，它们的值代表太阳的亮度。第一个数组是早上的测量值，第二个数组是晚上的测量值。实际上，我有大约80个这样的数组。我打算用matplotlib来绘制图像。绘制的圆圈在两种情况下大小是一样的，但由于地球的运动，图像的位置会稍微改变，这种变化是我想避免的。

>>> array1
[0, 0, 0, 0, 0, 0]
[0, 0, 0, 0, 0, 0]
[0, 0, 1, 3, 1, 0]
[0, 0, 1, 1, 2, 0]
[0, 0, 1, 1, 1, 0]
[0, 0, 0, 0, 0, 0]

>>> array2
[0, 0, 0, 0, 0, 0]
[0, 0, 0, 0, 0, 0]
[0, 0, 0, 0, 0, 0]
[0, 0, 1, 2, 1, 0]
[0, 0, 1, 1, 4, 0]
[0, 0, 1, 1, 1, 0]

在上面的例子中，较大的值表示亮的地方，而零值则被绘制成黑色的空白。两个数组的大小总是相同的。我该如何将第二个数组中重要的值（非零值）与第一个数组中的值对齐呢？所以最终的效果应该是这样的。

>>> array2(aligned)
[0, 0, 0, 0, 0, 0]
[0, 0, 0, 0, 0, 0]
[0, 0, 1, 2, 1, 0]
[0, 0, 1, 1, 4, 0]
[0, 0, 1, 1, 1, 0]
[0, 0, 0, 0, 0, 0]

这样做是为了能够有意义地对数组进行后续处理，比如计算平均值或总和等。需要注意的是，找到一个质心并相应对齐的方法不太有效，因为边缘可能会有高值，这些值在一天中会变化。

Answer 1

这种数据可能会出现的问题之一是，图像和像素的对齐并不完美。我用两个数组和其中的一个方块来说明我的观点：

array1:
0 0 0 0 0
0 2 2 2 0
0 2 2 2 0
0 2 2 2 0
0 0 0 0 0

array2:
0 0 0 0 0
0 1 2 2 1
0 1 2 2 1
0 1 2 2 1
0 0 0 0 0

如你所见，分辨率有限是个挑战，因为图像移动了0.5个像素。

当然，计算这两个数组的重心是很简单的，可以看到第一个数组的重心是(行, 列=2,2)，而第二个数组是(2, 2.5)。但如果我们把第二个数组向左移动0.5，就变成了：

array2_shifted:
  0   0   0   0   0
0.5 1.5 2.0 1.5 0.5
0.5 1.5 2.0 1.5 0.5
0.5 1.5 2.0 1.5 0.5
  0   0   0   0   0

这样一来，图像就开始分散开了。

当然，如果你的数组足够大，可能就不需要担心子像素的问题，但如果每个方向只有几或几十个像素，这可能会变得很麻烦。

解决这个问题的一种方法是先通过合适的外推来增大图像的大小（就像图像处理程序所做的那样；cv2模块提供了很多这样的功能）。然后可以以单像素的精度将图像对齐，再缩小回去。

---

无论如何，你需要一种方法来找出图像之间最佳的对齐位置。有很多选择需要做。一个重要的事情是，你可能不想把图像和第一个图像对齐，而是想把所有图像与一个参考对齐。在这种情况下，参考可以是图像中心的一个完美圆圈。这样你只需要移动所有图像以匹配这个参考。

一旦你选择了参考，就需要选择一种方法来获取一些关于图像对齐的指标。有几种可能性，但你可以从以下几种开始：

1. 计算图像的重心。

2. 计算图像与参考之间的相关性。相关数组中的最高点给你最佳匹配。

3. 在进行一些处理后再执行上述任一操作（通常是在每一端或两端限制动态范围）。

我会从这样的步骤开始：

• 可能需要对图像进行上采样（如果分辨率低的话）
• 限制动态范围的高端（例如 clipped=np.clip(image,0,max_intensity)）
• 计算重心（例如 scipy.ndimage.center_of_mass(clipped)）
• 根据重心的偏移量平移图像

对二维数组进行平移需要一些代码，但不应该太难。如果你确定周围是黑色的，可以使用：

translated = np.roll(np.roll(original, deltar, axis=0), deltac, axis=1)

这会将最左边的像素滚动到右边（或反之）。如果这样不好，那么你需要将它们清零。（或者看看：python numpy roll with padding）。

关于对齐过程的一个警告：最简单的方法（重心、相关性）在图像中有强度梯度时会失败。因此，你可能需要寻找边缘，然后再进行相关性计算。如果背景真的很黑，限制强度也会有所帮助。