Method ReadMidpointCoordinates

IO/ParallelNetCDF/vtkPSLACReader.cxx:947–1119

------------------------------------------------------------------------------

//------------------------------------------------------------------------------
int vtkPSLACReader::ReadMidpointCoordinates(
  int meshFD, vtkMultiBlockDataSet* vtkNotUsed(output), vtkSLACReader::MidpointCoordinateMap& map)
{
  // Get the number of midpoints.
  int midpointsVar;
  CALL_NETCDF_INT(nc_inq_varid(meshFD, "surface_midpoint", &midpointsVar));
  this->NumberOfGlobalMidpoints = this->GetNumTuplesInVariable(meshFD, midpointsVar, 5);
  if (this->NumberOfGlobalMidpoints < 1)
    return 0;

  vtkIdType numMidpointsPerPiece = this->NumberOfGlobalMidpoints / this->NumberOfPieces + 1;
  vtkIdType startMidpoint = this->RequestedPiece * numMidpointsPerPiece;
  vtkIdType endMidpoint = startMidpoint + numMidpointsPerPiece;
  endMidpoint = std::min(endMidpoint, this->NumberOfGlobalMidpoints);

  size_t starts[2];
  size_t counts[2];

  starts[0] = startMidpoint;
  counts[0] = endMidpoint - startMidpoint;
  starts[1] = 0;
  counts[1] = 5;

  vtkNew<vtkDoubleArray> midpointData;
  midpointData->SetNumberOfComponents(static_cast<int>(counts[1]));
  midpointData->SetNumberOfTuples(static_cast<vtkIdType>(counts[0]));
  CALL_NETCDF_INT(
    nc_get_vars_double(meshFD, midpointsVar, starts, counts, nullptr, midpointData->GetPointer(0)));

  // Collect the midpoints we've read on the processes that originally read the
  // corresponding main points (the edge the midpoint is on).  These original
  // processes are aware of who requested hose original points.  Thus they can
  // redistribute the midpoints that correspond to those processes that
  // requested the original points.
  std::vector<midpointListsType> midpointsToDistribute(this->NumberOfPieces);

  int pointsPerProcess = this->NumberOfGlobalPoints / this->NumberOfPieces + 1;
  for (vtkIdType i = 0; i < midpointData->GetNumberOfTuples(); i++)
  {
    double* mp = midpointData->GetPointer(i * 5);

    midpointPositionType position;
    position.coord[0] = mp[2];
    position.coord[1] = mp[3];
    position.coord[2] = mp[4];

    midpointTopologyType topology;
    topology.minEdgePoint = static_cast<vtkIdType>(vtkMath::Min(mp[0], mp[1]));
    topology.maxEdgePoint = static_cast<vtkIdType>(vtkMath::Max(mp[0], mp[1]));
    topology.globalId = i + startMidpoint + this->NumberOfGlobalPoints;

    // find the processor the minimum edge point belongs to (by global id)
    vtkIdType process = topology.minEdgePoint / pointsPerProcess;

    // insert the midpoint's global point id into the data
    midpointsToDistribute[process].position.push_back(position);
    midpointsToDistribute[process].topology.push_back(topology);
  }