如何在ParaView中修改VTKCompositeDataArray的对角元素？

问题背景 #

在ParaView可编程过滤器中，需导出名为metrique的VTKCompositeDataArray类型字段（每个点对应3×3矩阵，共9个值），目标是将数组中满足阈值条件的对角元素替换为指定值，但遇到错误：'VTKCompositeDataArray' object does not support item assignment。原代码对NumPy数组有效，但直接操作VTK数组失败；尝试手动转换VTKCompositeDataArray为NumPy数组时，又出现AttributeError: 'VTKCompositeDataArray' object has no attribute 'GetNumberOfArrays'错误。

解决方案 #

VTKCompositeDataArray是只读类型，无法直接进行元素赋值操作，需先将其转换为可修改的NumPy数组，处理完成后再转回VTK数组。以下是修正后的完整实现：

完整修改脚本 #

import numpy as np
import vtk
from paraview.util import numpy_support

# 参数配置
replacement_value = 0.5
lower_threshold = 0.5

# 读取输入数据
input0 = inputs[0]
Velocity = input0.PointData["Velocity"]

# 计算速度梯度与标量U
gradU = gradient(Velocity)
U = Velocity[:, 0]

# 构造张量B并扩展至所有点
B = np.array([[[0.1, 0, 0], [0, 0.1, 0], [0, 0, 0.1]]])
B_tensorial = np.tile(B, (input0.GetNumberOfPoints(), 1, 1))

# 将VTKCompositeDataArray转换为NumPy数组（关键步骤）
# 获取扁平化的底层VTK数组，再转为NumPy数组
gradU_vtk_flat = gradU.GetFlatArray()
gradU_np = numpy_support.vtk_to_numpy(gradU_vtk_flat)
# 重塑为(N, 3, 3)的张量形状
gradU_np = gradU_np.reshape((input0.GetNumberOfPoints(), 3, 3))

# 计算metrique的NumPy版本
metrique_np = gradU_np * B_tensorial * U[:, np.newaxis, np.newaxis]

# 替换对角元素中满足阈值条件的值
diagonal_indices = np.diag_indices(3)
diag_values = metrique_np[:, diagonal_indices[0], diagonal_indices[1]]
# 逻辑：绝对值大于阈值时替换为指定值，否则保留原数值
diag_values = np.where(np.abs(diag_values) > lower_threshold, replacement_value, diag_values)
metrique_np[:, diagonal_indices[0], diagonal_indices[1]] = diag_values

# 将处理后的NumPy数组转回VTK数组
metrique_flat = metrique_np.flatten()
metrique_vtk = numpy_support.numpy_to_vtk(metrique_flat, deep=True, array_type=vtk.VTK_DOUBLE)
# 设置数组名称与分量数（3×3矩阵对应9个分量）
metrique_vtk.SetName("metrique")
metrique_vtk.SetNumberOfComponents(9)

# 添加到输出点数据
output.PointData.append(metrique_vtk)

关键修正点 #

• 使用官方转换工具：借助paraview.util.numpy_support实现VTK与NumPy数组的互转，避免手动转换的兼容性问题。
• 正确访问VTKCompositeDataArray：通过GetFlatArray()获取其底层的扁平化VTK数组，而非错误调用GetNumberOfArrays()。
• 动态获取点数：用input0.GetNumberOfPoints()替代硬编码的nb_points，提升代码通用性。
• 恢复VTK数组结构：转回VTK数组时，需设置正确的分量数（9个），确保ParaView能识别为3×3张量字段。

内容的提问来源于stack exchange，提问作者Martin7