How to visualize vector fields on ImageData with Paraview?

[lahwaacz]

The examples directory contains several examples that show how to visualize data in ADIOS2 bp files with VTK and Paraview. However, all examples write only scalar fields and no vector field.

I thought the generalization would be to specify NumberOfComponents="3" same as in the XML VTK file formats. A small test code I used is:

#include <string>
#include <vector>

#include <adios2.h>

int main(int argc, char *argv[]){
    std::vector<double> data;
    std::size_t sizeR = 64;

    for(int x=0;x<sizeR;x++){
        for(int y=0;y<sizeR;y++){
            for(int z=0;z<sizeR;z++){
                data.push_back(10);
                data.push_back(20);
                data.push_back(30);
            }
        }
    }

    const std::string extent = "0 " + std::to_string(sizeR) + " " + "0 " + std::to_string(sizeR) + " " + "0 " + std::to_string(sizeR);

    const std::string imageData = R"(
            <?xml version="1.0"?>
            <VTKFile type="ImageData" version="0.1" byte_order="LittleEndian">
             <ImageData WholeExtent=")" + extent + R"(" Origin="0 0 0" Spacing="1 1 1">
                <Piece Extent=")" + extent + R"(">
                  <CellData>
                    <DataArray Name="vector_test" NumberOfComponents="3" />
                  </CellData>
                </Piece>
              </ImageData>
            </VTKFile>)";

    adios2::ADIOS adios;
    adios2::IO bpIO = adios.DeclareIO("bpIO");
    bpIO.SetEngine("BP4");

    adios2::Variable<double> var = bpIO.DefineVariable<double>("vector_test", {sizeR,sizeR,sizeR,3}, {0,0,0,0}, {sizeR,sizeR,sizeR,3});

    std::string filename = "vector-field-test.bp";
    adios2::Engine engine = bpIO.Open(filename, adios2::Mode::Write);

    engine.BeginStep();
    engine.Put(var, data.data());
    bpIO.DefineAttribute<std::string>("vtk.xml", imageData);
    engine.EndStep();
    engine.Close();

    return 0;
}

However, the vtkADIOS2VTXReader then complains:

terminate called after throwing an instance of 'std::runtime_error'
  what():  ERROR: NumberOfComponents 3 and variable names found 0 inside DataArray node Name in ADIOS2 VTK XML schema

It seems that the reader expects additional specification of components inside the DataArray element, but it accepts only pugi::node_pcdata nodes, i.e. plain-character nodes and no other elements. How is this supposed to be specified in the vtk.xml attribute?

Is there a working BP file with vector fields in an ImageData dataset that we could use to reverse-engineer the code for our writer? I found only this code in MFEM which writes correct vector fields, but for an UnstructuredGrid dataset.

[williamfgc]

@lahwaacz thanks for the detailed summary. I will try to find some time and take a look. Overall, you're correct. I remember implementing this for unstructured mesh, but I don't recall a use case for image data.

[lahwaacz]

@williamfgc Thanks, we would really appreciate it if you could take a look and give us some advice 🤞

[kmorel]

The current XML reader does not support vector data (as of ParaView 5.12.0). I think the best you can do right now is write out each component in independent fields and then use the Merge Vector Components filter to merge them together.

An alternate approach is to use the Fides reader, which is a newer schema based on JSON. The JSON schema is documented here: https://fides.readthedocs.io/en/latest/schema/schema.html.

Fides can read vector data directly. It is currently limited to sizing the mesh based on a scalar field, so in your example a dummy array of scalars sized to the regular mesh needs to be written. Here is a modification of your code to write out the data and an accompanying JSON file.

#include <fstream>
#include <string>
#include <vector>

#include <adios2.h>

int main(int argc, char *argv[]){
    std::vector<double> data;
    std::size_t sizeR = 64;

    for(int x=0; x<sizeR; x++){
        for(int y=0; y<sizeR; y++){
            for(int z=0; z<sizeR; z++){
                data.push_back(10);
                data.push_back(20);
                data.push_back(30);
            }
        }
    }

    std::vector<int> pointIndices;
    int index = 0;

    for(int x=0; x<sizeR+1; x++){
        for(int y=0; y<sizeR+1; y++){
            for(int z=0; z<sizeR+1; z++){
                pointIndices.push_back(index);
                ++index;
            }
        }
    }

    adios2::ADIOS adios;
    adios2::IO bpIO = adios.DeclareIO("bpIO");
    bpIO.SetEngine("BP4");

    adios2::Variable<double> var = bpIO.DefineVariable<double>("vector_test", {sizeR,sizeR,sizeR,3}, {0,0,0,0}, {sizeR,sizeR,sizeR,3});
    adios2::Variable<int> pointIdx = bpIO.DefineVariable<int>("point_indices", {sizeR+1,sizeR+1,sizeR+1}, {0,0,0}, {sizeR+1,sizeR+1,sizeR+1});

    std::string filename = "vector-field-test.bp";
    adios2::Engine engine = bpIO.Open(filename, adios2::Mode::Write);

    engine.BeginStep();
    engine.Put(var, data.data());
    engine.Put(pointIdx, pointIndices.data());
    engine.EndStep();
    engine.Close();

    std::ofstream jsonfile("vector-field-test.json");
    jsonfile << R"(
{
  "vector-field-test-data": {
    "data_sources": [
      {
        "name": "source",
        "filename_mode": "input"
      }
    ],
    "coordinate_system": {
      "array": {
        "array_type": "uniform_point_coordinates",
        "dimensions": {
          "source": "variable_dimensions",
          "data_source": "source",
          "variable": "point_indices"
        },
        "origin": {
          "source": "array",
          "values": [0, 0, 0]
        },
        "spacing": {
          "source": "array",
          "values": [1, 1, 1]
        }
      }
    },
    "cell_set": {
      "cell_set_type": "structured",
      "dimensions": {
        "source": "variable_dimensions",
        "data_source": "source",
        "variable": "vector_test"
      }
    },
    "fields": [
      {
        "name": "vector_test",
        "association": "cell_set",
        "array": {
          "array_type": "basic",
          "data_source": "source",
          "variable": "vector_test",
          "is_vector": "true"
        }
      }
    ],
    "step_information": {
      "data_source": "source"
    }
  }
})";
    jsonfile.close();

    return 0;
}

Here are the files it generates: vector-field-test.tar.gz

[lahwaacz]

@kmorel Thanks for the hint, I finally got to testing the Fides readers in Paraview. Unfortunately, it also does not fully meet all my expectations... ☹️

I've reproduced my issues with the Gray-Scott example and submitted ornladios/ADIOS2-Examples#90. One issue is that I apparently need at least one PointData in each dataset, which can be a waste of space when only CellData are used for the actual data. At some point, it would be even more efficient to use the rectilinear grid model instead of uniform image data... Another issue is the visual "gaps" between subdomains, any idea how to solve this?