def METHOD_NAME(heat, kernel=3):
    def max_pool2d(A, kernel_size, padding=1, stride=1):
        A = np.pad(A, padding, mode='constant')
        output_shape = ((A.shape[0] - kernel_size)//stride + 1,
                        (A.shape[1] - kernel_size)//stride + 1)
        kernel_size = (kernel_size, kernel_size)
        A_w = as_strided(A, shape=output_shape + kernel_size,
                         strides=(stride*A.strides[0],
                                  stride*A.strides[1]) + A.strides)
        A_w = A_w.reshape(-1, *kernel_size)
        return A_w.max(axis=(1, 2)).reshape(output_shape)
    pad = (kernel - 1) // 2
    hmax = np.array([max_pool2d(channel, kernel, pad) for channel in heat])
    keep = (hmax == heat)
    return heat * keep